A Method for Assessing Urban Ecological Resilience and Identifying Its Critical Distance Belt Based on the “Source-Sink” Theory: A Case Study of Beijing

Abstract

A reasonable assessment of urban ecological resilience (UER), as well as quantitative identification of critical thresholds of UER, is an important theoretical basis for the formulation of scientific urban development planning. The existing UER assessment methods ignore the dynamic relationship between protection factors and disturbance factors in urban systems and do not address the question of where UER starts to become unstable. Therefore, based on the “source-sink” landscape theory, we constructed a UER assessment model and a method to quantitatively identify the UER’s critical distance belt (UER-CDB) using the transect gradient analysis. Additionally, we combined scenario simulation to analyze the change characteristics of UER and its critical distance belt in different urban development directions over past and future periods. The results show that: (1) Based on the “source-sink” theory and transect gradient method, the UER can be effectively assessed and the UER-CDB can be quantitatively identified. (2) The UER in Beijing shows a distribution pattern of high in the northwest and low in the southeast, and the High resilience area accounts for more than 40%. (3) The changes in UER-CDB in Beijing in different development directions have obvious variability, which is mainly influenced by topography and policy planning. (4) Compared with the natural development scenario (NDS), the ecological protection scenario (EPS) is more consistent with Beijing’s future urban development plan and more conducive to achieving sustainable development. The methodology of this paper provides a fresh perspective for the study of urban ecological resilience and the critical threshold of ecosystems.

1. Introduction

Global urbanization is a significant trend in the development of human society in the 21st century [1,2], and the World Cities Report 2022: Envisaging the Future of Cities, published by UN-Habitat, indicated that 2.2 billion people will be added to urban areas in the next 30 years [3]. The vast population pressure will bring many ecological risks to cities, so building sustainable cities has become the most critical issue for human society in the 21st century [4,5]. This issue is even more critical in China, the world’s largest developing country. The urbanization rate has proliferated in the first two decades of the 21st century and exceeded the world average, providing new impetus for China’s sustained economic growth. However, blind urban expansion in the past has led to a series of ecological problems that threaten urban safety, such as the drastic reduction in natural resources [6], ecological fragility [7], environmental pollution [8], and the ability of cities to resist ecological risks has been decreasing, which has become an essential factor limiting sustainable urban development. Currently, in the context of China’s vigorous promotion of ecological civilization, the construction of healthier, safer, and more livable high-quality living spaces has become the primary goal of China’s urban development [9]. Therefore, how to reasonably assess ecological resilience, especially how to quantitatively identify the critical threshold of ecological resilience in urban systems, minimize the impact and loss of ecological risks to human society, and promote urbanization and natural ecosystems toward harmony and symbiosis, has become an important scientific problem that needs to be solved urgently [10].

The concept of resilience was first applied to ecology by the Canadian scholar Holling to characterize the ability of natural ecosystems to resist disturbance by human activities, restore equilibrium, and adapt to new environments [11]. With the increasing problems brought about by urbanization, the idea of resilience has been widely used by city managers and experts, and the concepts of “urban resilience” and a “resilient city” have been proposed one after another [12,13]. In 2016, the Third United Nations Conference on Housing and Sustainable Urban Development released the New Urban Agenda, which identified “resilient cities” as the core goal of future city construction [14]. Now, “urban resilience” and “resilient cities” have become a policy tool and strategic path for major cities worldwide to prevent and resolve significant risks. It has also become a research hotspot in sustainable urban development. Urban resilience refers to the ability of cities to withstand shocks, respond quickly, readjust, and recover more quickly in the face of disasters [15,16], and a resilient city is a city that has these capabilities.

Urban resilience is generally considered to be a highly complex coupled system consisting of an urban economy, society, institutions, ecology, and infrastructure [17]. The components are interdependent and promote each other. Among them, urban ecological resilience (UER) reflects the degree of coordination between urban development and ecosystems, and has a vital role in maintaining the homeostasis of urban ecological safety patterns [18], as well as serving as an external environmental safeguard for other resilient subsystems in the urban resilience system. In previous studies, many scholars focused on assessing the overall resilience level of cities [19]. In recent years, as the concept of urban ecological civilization has been raised, UER has become an essential tool for assessing the sustainable development of society, with the number of studies conducted for its separate evaluation and analysis gradually increasing. For example, Chen et al. [20] used urban ecological spatial land as a proxy variable for UER and argued that excessive construction land expansion in the urbanization process had exacerbated the UER loss in Guangzhou. Some scholars also tried quantitatively assessing UER based on geospatial analysis methods. For instance, Wang et al. [21] constructed an assessment system of UER based on the essential characteristics of ecological resilience with three subsystems of scale, density, and morphology. In conclusion, most of the existing studies have assessed UER based on the characteristics of the ecosystem itself.

It has been argued that UER results from the interaction between human needs and the natural environment. Its magnitude is determined by various protection and disturbance factors inside and outside the urban system [10,22]. Previous studies on UER have only considered the protection factors in ecosystems and neglected the driving role of disturbance factors on ecological resilience. However, the “source-sink” theory in landscape ecology provides an idea to solve the above problems. According to the concept proposed by Chen et al. [23], “source” refers to the factors that promote the development of ecological processes, and “sink” refers to the factors that prevent and delay the development of ecological processes, which correspond to the protection and disturbance factors in urban ecosystems, respectively. According to the self-organization criticality theory [24] and the adaptive cycle theory of resilience [25], the magnitudes of “source” and “sink” factors in urban ecosystems change continuously with the natural socio-economic development, and there are prominent characteristics of urban–rural gradient changes. At the periphery of the built-up area, the “source” is stronger than the “sink”, and the urban ecosystem is relatively stable. When close to the built-up area, where the difference between the two forces gradually decreases and the system resilience approaches a critical threshold, the ecosystem accelerates from stable to unstable [26,27]. The assessment of ecological resilience has been initially discussed by scholars based on the “source-sink” theory [22,28]. However, the critical threshold of ecological resilience has not been further explored, which means that the question of where ecological resilience starts to move from a stable to an unstable state is not addressed. According to Holling [11], the ecological threshold includes threshold points and threshold belts. In normal urban development, the “source” factor and “sink” factor are a slow interaction process which is more in line with the characteristics of ecological threshold belts. Therefore, based on the above theoretical basis and the actual development of cities along the urban–rural gradient, this paper concretizes the concept of the critical threshold of UER: The critical threshold of UER can be expressed as the maximum distance from the area where the “source” and “sink” factors reach equilibrium to the urban center while maintaining ecosystem homeostasis [20]. We refer to it as the urban ecological resilience critical distance belt (UER-CDB).

In 2021, “resilient cities” were officially included in China’s 14th Five-Year Plan and 2035 Vision [6]. Therefore, the next decade will be an important period for China to construct “resilient cities” and realize sustainable urban development. The main objective of this paper is to explore a method to assess the UER and quantitatively identify the UER-CDB based on the “source-sink” theory, in order to intuitively reflect the dynamic relationship between urbanization and ecosystems in the past 20 years and under the natural development scenario (NDS) and ecological protection scenario (EPS) in 2030. This paper takes Beijing, the capital city of China, as a case study, and the specific research contents are as follows:

• Assessment of the spatial and temporal characteristics of land use change from 2000 to 2020 and under different scenarios in 2030.
• Assessment of the spatial and temporal characteristics of UER change from 2000 to 2020 and under different scenarios in 2030.
• Assessment of the change characteristics of UER-CDB from 2000 to 2020 and under different scenarios in 2030.

This paper not only provides an important scientific reference for urban land spatial planning and resilient city construction, but also the method of assessment of UER and quantitative identification of UER-CDB using the “source-sink” theory and its application to the historical period and different development scenarios in the future, which will be a new attempt in the study of the ecological resilience and critical threshold of ecosystems.

2. Materials and Methods

2.1. Study Area

Beijing is located in the northern part of the North China Plain (115°25′–117°30′E, 39°28′–41°05′N), with a total area of about 16,410 square kilometers. It is adjacent to Tianjin City and Hebei Province, forming the Beijing–Tianjin–Hebei region together. The topography of Beijing is high in the northwest and low in the southeast. The western and northern parts are surrounded by mountains, with elevations between 1000 and 1500 m, and are mainly covered by woodlands and grasslands. The plain area in the southeast is about 6338 km² and is the main concentration of urban development. Beijing has a warm, temperate, semi-humid, and semi-arid monsoon climate with an average annual temperature of 12.3 °C and an annual precipitation of about 600 mm. Beijing is the political center, cultural center, international communication center, science and technology innovation center of China, and a megacity, with 16 districts under its jurisdiction: Dongcheng, Xicheng, Haidian, Chaoyang, Fengtai, Shijingshan, Mentougou, Fangshan, Tongzhou, Shunyi, Changping, Daxing, Huairou, Pinggu, Yanqing, and Miyun (Figure 1).

The first 20 years of the 21st century were crucial for Beijing to promote urbanization. Beijing’s resident population in 2020 reached 21.893 million, an increase of 60% compared with 2000, and its GDP in 2020 was CNY 3610.26 billion, 14 times higher than that in 2000. Under the combined influence of the natural environment and social economy, Beijing’s urban land has gradually expanded over the past 20 years, with a development intensity of about 45% in plain areas. The conflict between ecological security and urban development has become increasingly prominent, and local governments are faced with the challenge of balancing ecological protection and urban development [29]. To alleviate this conflict, the Beijing Municipal Government has formulated several ecological conservation-related plans to strictly control the increase in construction land and gradually transition from stock development to reduced development. Among them, the Special Plan for Beijing’s Ecological Security Pattern (2021–2035) clearly states that it is necessary to build a resilient urban ecosystem, effectively guarantee the ecological security of the capital, improve the ecological quality of the city, and continuously meet the residents’ demand for high-quality ecological space. Ultimately, Beijing will be built into an international eco-city model. Based on this, this paper took Beijing as the study area to explore the dynamic changes in UER and UER-CDB during the urbanization process, which has representative significance.

2.2. Dataset and Preprocessing

We took full account of the availability and applicability of the data and used the following data:

(1) Land use data: we compared three data products, including Globeland30, CLCD, and GLC_FCS30, and finally chose Globeland30 (http://globeland30.org/, accessed on 30 July 2022) as the land use data source because the time scale and data quality of these data fit our research needs. First, this dataset contains three times—2000, 2010, and 2020—and effectively reflects the global land use changes over the past 20 years, which is valuable for monitoring global environmental changes and regional resource management [30,31]. Second, these data have a spatial resolution of 30 m and are generated based on satellite images such as Landsat-5, Landsat-8, HJ-1 (China Environment and Disaster Mitigation Satellite), and GF-1 (China High Resolution Satellite) [32]. The total accuracy of the GlobeLand30 V2010 and V2020 is 83.50% and 85.72%, respectively, with Kappa coefficients of 0.78 and 0.82. These results were obtained from the validation of more than 150,000 points and 230,000 points, respectively [29]. Additionally, we found that the data better distinguished between cultivated and forested land than the other two data, and the patches were more aggregated. In addition, to make these data more realistic, we combined them with the manual correction performed by Google Images in Beijing.

We reclassified these data into six categories: cropland, forest, grass, water, building, and bare land, as the basic data for land use change analysis and simulation of the spatial distribution of land use in 2030.

(2) Land use change driving factor data: These data mainly include two aspects, socioeconomic and natural climate, as the driver data for simulating the spatial distribution of land use in 2030 (

Table 1. Data on drivers of land use change.

Category	Data	Resolution	Data Resource
Socioeconomic driver data	GDP	1 km	https://www.resdc.cn, accessed on 1 August 2022
	Population density	100 m	https://hub.worldpop.org/, accessed on 1 August 2022
	Nighttime lights	500 m	https://eogdata.mines.edu/products/vnl/, accessed on 1 August 2022
	Road network		Ministry of Natural Resources of China
	Governments point		Ministry of Natural Resources of China
Nature and climatic driver	DEM	30 m	ASTER GDEM 30 M dataset http://www.gscloud.cn, accessed on 1 August 2022
	Temperature	1 km	http://www.geodata.cn, accessed on 1 August 2022
	Precipitation	1 km	http://www.geodata.cn, accessed on 1 August 2022
	River system	2 m	Ministry of Natural Resources of China

). The socioeconomic data include GDP grid data, WorldPop global population density data, NPP/VIIRS annual nighttime light data, and data on various levels of roads and government points from the results of China’s basic geographic national monitoring. Based on the above data, ten factors such as nighttime lights, population density, GDP, distance from national roads, distance from provincial roads, distance from urban roads, distance from county roads, distance from railroads, distance from railroad stations, and distance from county government points were selected and calculated as the socioeconomic drivers of land use change. Natural climate data include DEM, temperature, precipitation, and river system data from China’s basic geographic national monitoring results. Based on the above data, five factors, such as elevation, slope, temperature, precipitation, and distance from water bodies, were selected and calculated as natural climate drivers of land use change.

(3) Government planning data: Beijing Urban Master Plan (2016–2035); Beijing 14th Five-Year Plan for Land Resources Protection and Utilization; Beijing Special Plan for Ecological Security Pattern (2021–2035), used for extracting Beijing’s ecological protection red line boundary, cropland protection red line boundary, and other information, were used as the restriction data for the simulation of land use in 2030.

The projection coordinates of the above data were uniformly transformed to CGCS2000_3_Degree_GK_CM_117E, and the resolution of the raster data was uniformly transformed to 30 m × 30 m.

2.3. Technological Process

The detailed technical process of this paper includes the following steps (Figure 2):

(1) Simulating the spatial distribution of land use in 2030. Based on the Markov matrix and the patch-generating land use simulation (PLUS) model, we simulated the spatial distribution of land use in 2030 under the NDS and EPS and analyzed their spatial and temporal changes.

(2) Assessing the UER based on the “source-sink” theory. Firstly, resilience was divided into two parts: basic resilience and process resilience. Then, the resilience values of the “source” and “sink” factors were calculated based on the “source-sink” landscape theory, and then the UER was calculated comprehensively. Finally, we analyzed the spatial and temporal changes in UER.

(3) Identifying the UER-CDB. Firstly, eight transects were set up in eight directions: east, northeast, north, northwest, west, southwest, south, and southeast. Then, the “source” resilience values and “sink” resilience values on the transects were extracted according to the distance gradient and fitted into curves. Finally, the critical UER-CDB was identified according to the intersection of the two curves.

2.4. Simulating the Spatial Distribution of Land Use in 2030

It has been shown that UER can be effectively improved through rational urban planning [33]. Additionally, the prerequisite of planning is to predict the change characteristics of urban land use and its ecological effects under different development strategies, and simulating the development characteristics of future land use based on different scenarios is an important method to achieve the above goal [34]. In this paper, based on the characteristics of land use changes in Beijing in historical years and its future regional spatial development plan, we designed two development scenarios (NDS and EPS). Then, we used the PLUS model and Markov model to simulate the future land use distribution, and then assessed the future UER and identified UER-CDB on this basis.

2.4.1. Setting Development Scenarios and Predicting Land Use Demand

(1)

NDS

Based on the probability of land use transfer in the study area from 2010 to 2020, we used the Markov model to predict the land use demand in 2030 and did not set any restriction zone in the simulation process. The Markov model is a pixel-scale-based model that generates a land use transfer matrix from land use data in different periods and calculates the transfer probability of land use change in the study area matrix, which has a robust quantitative prediction capability [35] and has been widely used for land demand prediction [36]. The specific process is expressed as follows:

$$S_{t+1}=P_{ij}\times S_i$$

(1)

where S_t and S_t+1 denote the land use status at t and t + 1, respectively; P_ij is the land use transfer probability matrix, which denotes the probability of transferring land type i to land type j.

(2)

EPS

Based on the land use transfer probability under the NDS and with reference to the Beijing 14th Five-Year Plan for Land Resources Protection and Utilization planning objectives, the conversion probability of ecological space land to the building was strictly limited. The scale of the building was controlled within 3670 km², while the ecological protection red line and the cropland protection red line were used as restriction zones in the simulation process.

2.4.2. The PLUS Model and Accuracy Verification

The PLUS model was constructed based on the rule mining framework of the Land Expansion Analysis Strategy (LEAS) and the CA model of the multi-type random patch seeds (CARS) [37]. Among them, LEAS focuses on analyzing the degree of contribution of different drivers to each land use type change one by one through random forest classification (RFC) and then determining the probability of future expansion of different land use types in the study area [38]. In this paper, we selected a total of 15 drivers, both socio-economic and natural climate (see “Section 2.2”). The CARS model is based on the probability of expansion of each land use type, with constraints of adaptive coefficients and neighborhood effects, to drive land use to meet future demand [39]. Compared with other models, the PLUS model simplifies the analysis of land use change while maintaining the ability to support multiple types and complex land use changes. Moreover, it can better simulate changes at the patch level of natural ecological land, such as forests and grasses, with higher accuracy [40].

To verify the accuracy of the PLUS model in simulating the land use distribution in Beijing, we first simulated the land use in 2020 based on the land use data in 2000 and 2010 and 15 drivers. Then, we compared the simulation results with the actual land use distribution in 2020 by calculating the Kappa and the figure of merit (FoM) to verify the simulation accuracy. Among them, the Kappa is calculated as follows:

$$K=\frac{P_a-P_b}{P_c-P_b}$$

(2)

where K represents Kappa value, P_a is the proportion of correct simulations, P_b is the proportion of correct model simulations in the stochastic case, and P_c is the proportion of correct simulations in the ideal case, usually defined as 1. The Kappa takes values in the range [−1, 1], and it is usually considered that if the Kappa is higher than 0.75, the model simulation results achieve a high level of agreement with the actual distribution [41]. The FoM is more advantageous in measuring the goodness-of-fit in land use change simulations [42,43]. It is calculated as follows:

$$F=\frac{1}{\max \left\{N_e,N_d\right\}}{\displaystyle \sum _{k=1}^{N_{\mathrm{d}}}\frac{1}{\left(1+\beta d{\left(k\right)}^2\right)}}$$

(3)

where F represents FoM value, N_e denotes the pixel volume of the land use simulation result; N_d denotes the actual pixel volume of the land use; β is a scale factor greater than 0, usually 1/9; and d(k) is the distance between the kth detected pixel of the actual land use and the pixel of the simulation result. Generally, the FoM value is within 0.3, while in practice, FoM usually takes values between 0.1 and 0.2, indicating a relatively high accuracy [44]. The validation results show that when the sampling rate is 5%, the Kappa is 0.753, and the FoM is 0.181, indicating that the PLUS model simulates land use in Beijing with high accuracy suitable for the subsequent study.

2.5. Assessing UER Based on “Source-Sink” Theory

Referring to the theory proposed by Wang et al. [10], resilience can be divided into basic resilience and process resilience. Among them, basic resilience characterizes the level of resilience in the basic configuration state of the system, and process resilience characterizes the level of resilience when the system is subjected to external disturbances. In addition, to fully reflect the spatial heterogeneity and continuity of resilience, we set a 1 km cell grid as the primary research unit based on the minimum patch area of available land use data. We calculated the resilience index in these units.

2.5.1. Basic Resilience

Basic resilience can be divided into three levels of resilience in terms of scale, morphology, and layout [21,45]. Based on the theory of landscape ecology and the definition of ecological space in the Beijing Special Plan for Ecological Security Pattern (2021–2035), we considered cropland, forest, grass, water, and bare land in land use data as “source”. Then, we regarded the building as “sink”. We calculated the area percentage, landscape shape index, and nearest distance index of the “source” and “sink” factors in the 1 km grid cell, respectively. Finally, the indices were weighted using the entropy method and linearly summed to obtain the basic resilience of the “source” and “sink” factors, respectively. The specific calculation formula is as follows:

$${\mathit{BR}}_{\mathit{source}}({\mathit{BR}}_{\mathit{sink}})=\mathit{AP}\times 0.328+\mathit{LSI}\times 0.203+\mathit{NDI}\times 0.469$$

(4)

where BR_source and BR_sink are the basic resilience of the “source” and “sink” factors, respectively, and AP, LSI, and NDI are the normalized values of the area percentage, shape index, and nearest distance index of the “source” and “sink” factors in each cell, respectively. Since the data of these three indices are discrete and have no primary and secondary relationships with each other, the average weights of the three indices calculated using the entropy weight method are 0.328, 0.203, and 0.469, respectively.

2.5.2. Process Resilience

The ability of ecosystems to resist external disturbances has been shown to be related to ecosystem service functions [33,46,47]. Ecosystem service supply and demand together form a dynamic process of ecosystem service flow between urban systems, and the combination of ecosystem service supply and demand assessment can accurately reflect the interaction of protection and disturbance factors in the system [48]. Therefore, in this study, based on the relationship between the supply and demand of ecosystem services, the matrix method [49] was used to assign values to each land use type according to the actual development of Beijing. Additionally, the literature results [50] were referred to in order to obtain the supply and demand matrix of ecosystem services by assigning values to each land use type. Then, the supply side of ecosystem services was set as “source” with the demand side was set as “sink,” and the process resilience of “source” and “sink” factors were measured separately. The calculation formula is as follows:

$${\mathit{PR}}_{\mathit{source}}({\mathit{PR}}_{\mathit{sink}})={\displaystyle \sum _{i=1}^{\mathrm{n}}{\displaystyle \sum _{j=1}^{\mathrm{m}}\left(A_i\times P_{ij}\right)}}$$

(5)

where PR_source and PR_sink are the process resilience of “source” and “sink” factors, respectively, A_i is the area of the ith land use type, and P_ij is the score of supply and demand of the jth ecosystem service of the ith land use type.

2.5.3. Ecological Resilience

The urban system is a dynamically changing open system, and basic resilience and process resilience play an equally important role in the system’s resilience. Therefore, in this study, the basic resilience and process resilience weights were set to 0.5, and the “source” resilience and “sink” resilience were calculated using linear summation. The calculation formula is as follows:

$$R_{\mathit{source}}(R_{\mathit{sink}})={\mathit{BR}}_{\mathit{source}}({\mathit{BR}}_{\mathit{sink}})\times 0.5+{\mathit{PR}}_{\mathit{source}}({\mathit{PR}}_{\mathit{sink}})\times 0.5$$

(6)

where R_source is the “source” resilience and R_sink is the “sink” resilience. To facilitate comparative analysis, this paper further drew on the concept of resilience in economics and used the ratio of the difference between the two and the “sink” resilience as the UER, as follows:

$$\mathit{UER}=\frac{R_{\mathit{source}}-R_{\mathit{sink}}}{R_{\mathit{sink}}}$$

(7)

Resilience has a certain temporal effect, and although the concept of resilience in this paper expressed a time-slice result, we tried to analyze the temporal and spatial dynamics of ecological resilience in the process of urban development by calculating its values at key time points. Therefore, to facilitate the statistical analysis of its dynamic change pattern, we divided the UER into five levels. When UER < 0, the “source” resilience is smaller than the “sink” resilience, and the UER level is the lowest, so it is defined as a Low grade. When the UER is ≥0, it is divided into four levels according to the equating method, namely Lower (0–0.25), Moderate (0.25–0.5), Higher (0.5–0.75), and High (0.75–1).

2.6. Identifying the UER-CDB

Since urban planning and development are directional, in order to reflect the regular changes in urban development and ecological patterns along specific directions in the study area, scholars have mainly used the transect gradient method to demonstrate the spatial distribution patterns of the study targets [51,52]. This method can reveal the spatial heterogeneity of urbanization and ecological patterns in the process of urban development and their formation mechanisms. Therefore, in this study, we analyzed the UER-CDB in different development directions in Beijing using the transect gradient method. We referred to previous studies [53] and set up transects along the main development axes of Beijing based on the Beijing Urban Master Plan (2016–2035). In this plan, there is a great similarity between the urban development direction and the geographical direction of Beijing. Therefore, transect directions can be divided into the following three types:

• The extension of the Beijing central axis as the north–south direction and the extension of Chang’an Street as the east–west direction, which are the significant development axes in Beijing’s urban planning.
• The intersection of Beijing’s central axis and Chang’an Street as the urban center, the extension of the connection line between the urban center and Miyun District Government as the northeast direction, and the extension of the connection line between the urban center and Yanqing District Government as the northwest direction. These two directions are located on the connection line between Beijing’s urban center and ecological cultured area, and are also essential areas for promoting the integrated development of urban and rural areas.
• The extension of the connection line between the urban center and Fangshan District Government as the southwest direction, and the extension of the connection line between the urban center and Yizhuang New Town (Yizhuang Town Government) as the southeast direction. These two directions are located on the “Beijing-Baoding-Shijiazhuang” development axis and the “Beijing-Tianjin” development axis, respectively, which are essential areas reflecting the coordinated development of the Beijing–Tianjin–Hebei region.

Therefore, the transects were finally set in eight directions: north, northeast, east, southeast, south, southwest, west, and northwest (Figure 3).

Studies have shown that transect gradient analysis tends to be scale-dependent, with transect bandwidth and step size being crucial factors [54]. Therefore, to fully reflect the gradient regularity of the “source” and “sink” resilience on the transect and to avoid interference caused by local differences, a suitable transect bandwidth and step size should be chosen for the gradient analysis. In setting the transect bandwidth, the “source” and “sink” factors were considered to have a spillover effect, which not only affected their environment but also had a radiative effect on the surrounding environment [22]. Since the old city of Beijing is approximately a square of 8 km×8 km, and the subsequent urban development is centered on the old city and expands outwards, the width of the transect was set at 8 km, which can fully reflect the influence and change pattern of “source” resilience and “sink” resilience in the process of urban development. In setting the transect step size, we have considered two aspects. First, according to the definition of UER-CDB in this paper, the distance should be a distance interval rather than a particular distance. Secondly, when the step size is small due to the interference of local features, the resilience index fluctuates sharply along the transect, and the changing pattern is not apparent. When the step size is large, it will cover some vital information and cannot correctly reflect the spatial variation of resilience on the transect [54]. After several experiments, we found that a step size of 5 km could fully reflect the spatial differentiation of resilience without causing the curve to fluctuate drastically.

Therefore, in this study, rectangular quadrats of 5 km × 8 km were finally set up, and the average values of the “source” and “sink” resilience in each quadrat were calculated. Then, the distance was used as the horizontal axis, and these two-factor resilience values were used as the vertical axis to plot the curves. Finally, we determined the distance interval where the intersection point of the two curves is located. Moreover, it should be specially noted that even if the quadrat size is set to 5 km × 8 km, the “source” and “sink” resilience do not vary strictly according to the gradient in some directions during the actual development process of the city. The relationship between the “source” and “sink” resilience values might change abruptly when fitting the curves, for many reasons, such as ecological parks within built-up areas and airports in suburban areas, or the development of suburban towns. This will result in multiple intersections of the two curves, so we needed to remove these intersections due to mutations and use the distance interval where the remaining intersections are located as UER-CDB (Figure 4).

3. Results

3.1. Analysis of Spatiotemporal Variation Characteristics of Land Use

From 2000 to 2020, the largest area of each land use type in Beijing was forest land, followed by cropland, building, and grass, with only a tiny proportion of water and bare land (Figure 5a–c). The spatial distribution of each land use type was strongly influenced by topography. The mountain area in the west and north was unsuitable for large-scale development. It was therefore dominated by forest cover, which increased slightly over the 20 years, from 43.67% to 45.08% (Figure 5f). The central and southeastern areas of the city were flat, with more cropland and building and frequent human activities. The years 2000–2020 saw a marked expansion of the building, with the area share increasing from 9.69% to 21.20%. The expansion was more significant in the last decade than in the first, with the expansion area gradually extending from the central districts of Haidian, Chaoyang, and Fengtai to the outer districts of Changping, Shunyi, Tongzhou, Fangshan, and Daxing. The expansion of the building has also encroached on a large area of cropland, resulting in a 10.77% reduction in the proportion of cropland. Beijing’s water was mainly located in the Miyun Reservoir, the Yongding River Basin, and the Chaobai River Basin, and its area has fluctuated over the past 20 years. The grass area was mainly located in the gently sloping mountainous areas, the mudflat areas around the waters. Additionally, the bare land area was relatively small and mainly distributed in the mountainous areas.

The spatial distribution of the various land use types under the two development scenarios of NDS and EPS in 2030 shows significant differences, with the most apparent changes still in the central city and its peripheral areas. In the NDS (Figure 5d), the building area will reach 4328.62 km², an increase of 849.82 km² compared with 2020 and much higher than the planned limit of 3670 km². In the south and southeast, the building area will occupy most of the plains of Tongzhou, Fangshan, and Daxing. Additionally, in the north, the building will gradually extend to Miyun, Huairou, and Pinggu through Changping and Shunyi. At the same time, the area of cropland will decrease significantly, by 887.65 km² compared with 2020, and the proportion of the total area will start to be lower than that of the building, which is only 18.39%. However, the trend of changes in other land uses will be more moderate. In the EPS (Figure 5e), the forest will increase by 125.49 km² and the water area by 77.08 km² compared with 2020. In addition, the trend of building expansion and cropland reduction will slow down. Compared with the NDS, the building area will decrease by 737.26 km² and remain within 3670 km², and cropland will increase by 653.83 km². It can be seen that under the NDS, as urbanization progresses, the urban building area will expand significantly and continue to encroach on a large area of basic agricultural land and ecological space, and the sustainability of urban land resources will face a considerable challenge. Under the EPS, however, the building area will slow down due to the strict protection of permanent basic agricultural land and ecological space, thus maintaining the sustainable development of urban land resources.

3.2. Analysis of Spatiotemporal Variation Characteristics of UER

From 2000 to 2020, the spatial distribution characteristics of UER were consistent with the spatial distribution pattern of land use. Among them, the High and Higher resilience were mainly distributed in the western and northern mountainous areas, the Moderate and Lower resilience were mainly distributed in the periphery of the central and southeastern plain areas, while the Low resilience was mainly located in the urban building areas (Figure 6a–c). In terms of temporal changes, the UER in Beijing changed significantly over the 20 years, with the changing area accounting for 34.54% of the total area and dominated by a decrease in level (Figure 6f,

Table 2. Transfer matrix for different levels of UER from 2000 to 2020.

Level	Low	Lower	Moderate	Higher	High	Total 2000
Low	1445.59	125.45	23.16	2.89	0.00	1597.09
Lower	1219.77	1203.37	414.96	45.36	2.89	2886.35
Moderate	447.76	358.98	689.99	41.49	6.76	1544.98
Higher	410.13	744.99	465.13	822.19	127.38	2569.82
High	32.81	138.00	788.41	271.17	6581.37	7811.76
Total 2020	3556.06	2570.79	2381.65	1183.10	6718.40	16,410.00

Note: Unit of area—km².

). Among them, the Low resilience area kept expanding outwards, with the area share increasing from 9.73% in 2000 to 21.67% in 2020, mainly by the transfer of Lower resilience, with 1219.77 km², followed by Moderate (447.76 km²) and Higher (410.13 km²). It showed a gradual trend from scattered polycentric to overall coalescence, in line with building area expansion. The High and Higher resilient areas also decreased from 63.26% in 2000 to 48.15% in 2020, with Higher resilience mainly transforming into Lower resilience (744.99 km²) and High resilience mainly transforming into Moderate resilience (788.41 km²). This situation has led to the gradual fragmentation of High resilient areas in the western and northern mountainous regions from being fully connected, suggesting that the rapid expansion of the building area has dramatically affected the level of UER of the surrounding areas.

Under both development scenarios in 2030, Beijing’s UER will remain dominated by the High and Higher classes, while the spatial variability of the other levels will be more pronounced (Figure 6d–f). Compared with 2020, the High and Higher resilience areas under the NDS have relatively little change, while Low resilience will expand by 38.86% to an area of 4937.96 km² in 2030, of which 1166.70 km² are from the transfer of Lower resilience, followed by that of Moderate resilience (

Table 3. Transfer matrix for different levels of UER under two scenarios from 2020 to 2030.

Level	Low	Lower	Moderate	Higher	High	Total 2020
Low	3485.61	68.52	1.93	0.00	0.00	3556.06
Lower	1166.70	1167.66	220.02	16.41	0.00	2570.79
Moderate	284.68	712.18	1225.56	132.21	27.02	2381.65
Higher	0.97	87.81	205.55	863.68	25.09	1183.10
High	0.00	0.97	64.65	39.56	6613.22	6718.40
Total 2030 (NDS)	4937.96	2037.14	1717.71	1051.86	6665.33	16,410.00
Low	3443.16	111.94	0.96	0.00	0.00	3556.06
Lower	441.02	1867.29	250.90	11.58	0.00	2570.79
Moderate	8.69	353.19	1943.53	70.45	5.79	2381.65
Higher	0.00	6.74	49.22	1072.13	55.01	1183.10
High	0.00	0.00	12.55	4.82	6701.03	6718.40
Total 2030 (EPS)	3892.85	2339.18	2257.16	1158.98	6761.83	16,410.00

Note: Unit of area—km²; NDS—Natural Development Scenario; EPS—Ecological Protection Scenario.

). Low resilience will occupy most of the plains and increase significantly in the northwest around the Yanqing District, while Lower and Moderate resilience areas will decrease significantly. Under the EPS, High and Higher resilience areas will show an increasing trend, while Low resilience areas will decrease significantly, and only 441.02 km² of Lower resilience and 8.69 km² of Moderate resilience transferred in. The spatial distribution pattern of Lower and Moderate resilience will be the same as that in 2020.

According to the natural development trend, the building area in Beijing will further expand in 2030, and the “sink” factors in the plains will tend to converge, resulting in a decrease in the connectivity and stability of the “source” factors in the surrounding areas, which will further lead to a significant increase in the Low resilience area. Due to topographic constraints, the High and Higher resilience areas remain stable in the western and northern mountainous areas. However, the fragmentation situation has yet to improve, and the ecological risk will remain serious. The EPS restrains the expansion of the building area, enhances the integration of the spatial structure of the new buildings, and slows down the increase in landscape fragmentation and separation in urban fringe areas due to the uncontrolled expansion of the building area so that the UER in Beijing in 2030 can remain stable.

3.3. Analysis of Variation Characteristics of UER-CDB

We extracted the “source” and “sink” resilience values in different directions and plotted the curves as shown in Figure 7 and Figure 8. In the north, east, southeast, and west, the “source” resilience tends to increase as the distance to the city center increases, while the “sink” resilience tends to decrease, with only one intersection point between the two curves. In other directions, the change in the “source” and “sink” resilience does not follow a strict urban gradient, so there are multiple points of intersection.

From 2000 to 2020, the UER-CDB in Beijing varied significantly in different development directions. The UER-CDB in the north gradually increased outward from 15–20 km to 30–35 km over 20 years (Figure 7a,

Table 4. The UER-CDB in eight directions from 2020 to 2030.

Direction	2000	2010	2020	2030 (NDS)	2030 (EPS)
North	15–20	20–25	30–35	35–40	30–35
Northeast	10–15	10–15	45–50	50–55	45–50
East	15–20	30–35	30–35	—	30–35
Southeast	10–15	20–25	25–30	35–40	25–30
South	5–10	20–25	20–25	30–35	20–25
Southwest	15–20	35–40	40–45	55–60	40–45
West	20–25	20–25	25–30	25–30	25–30
Northwest	15–20	20–25	45–50	45–50	45–50

Note: Unit—km; NDS—Natural Development Scenario; EPS—Ecological Protection Scenario.

). For the northeastern area of the city, the UER-CDB was 10–15 km from 2000 to 2010, while the development of the Capital International Airport caused a sudden change in the “source” and “sink” resilience of the surrounding area, resulting in other intersections of the two curves at around 25 km. In 2010, the impact of the airport was further strengthened and the development of the Miyun District led to new intersections at around 65 km. By 2020, the UER-CDB to the northeast increased to 45–50 km, and the former UER-CDB became an inner-city ecological park. Additionally, the impact of the Miyun District on the UER of the surrounding area continued to increase (Figure 7b, Table 4). The eastern topography is dominated by plains, so there was a significant expansion of the building area in the past 20 years, and the “source” resilience remained low. Additionally, the UER-CDB increased from 15–20 km in 2000 to 30–35 km in 2020 (Figure 7c, Table 4). In the southeastern direction, which was also a plain area, the UER-CDB expanded from 10–15 km in 2000 to 25–30 km in 2020 (Figure 7d, Table 4).

For the southern part of the city, the UER critical distance expanded from 5–10 km in 2000 to 20–25 km in 2010. Moreover, by 2020, although the critical distance was still 20–25 km, the construction and operation of Daxing International Airport affected the ecological resilience of the surrounding area, resulting in a new intersection of the two curves at around 45 km (Figure 7e, Table 4). In the southwestern direction, the UER-CDB was 15–20 km in 2000, and the development of the Fangshan District led to a sudden change in the “source” and “sink” resilience at 30 km. In 2010, the UER-CDB increased to 35–40 km and the building area in the Fangshan District was gradually connected to the central city. By 2020, the UER-CDB further increased to 40–45 km (Figure 7f, Table 4). In the west of the city, there is a clear demarcation between plain and mountainous areas, so the “source” and “sink” resilience varied significantly along the topographic gradient, with the UER-CDB increasing from 20–25 km to 25–30 km (Figure 7g, Table 4). For the northwest of the city, the plain area alternates with the mountainous area, so the curve fluctuated considerably. In 2000, the UER-CDB was 15–20 km, but at 40 km, the “source” and “sink” resilience changed abruptly due to the development of the Changping District. In 2010, the UER-CDB increased outwards to 15–20 km. Additionally, by 2020, the building area of the Changping District gradually connected with the central city, resulting in a further increase in the UER-CDB outwards to 45–50 km, reaching the northwestern mountains. At the same time, in Yanqing, 70 km to the northwest, the gap between the “source” and “sink” resilience gradually decreased from 2000 to 2020, and new intersections will soon be formed (Figure 7h, Table 4).

In the NDS of 2030, the UER-CDB in the west and northwest will remain unchanged compared to those of 2020 due to the constraints of the mountainous topography, while the UER-CDB in the other directions will increase further outwards and will be more significant in the direction dominated by the plains (Figure 8, Table 4). In particular, the eastern building area will be gradually connected to Hebei Province, so there will be no UER-CDB within Beijing by 2030. In the southwest, the UER-CDB will increase by around 15 km, followed by the southeast and south, which will increase by around 10 km and gradually approach Hebei Province, while the other directions will increase by around 5 km. In addition, the development and construction of the Miyun District in the northeast, Daxing International Airport in the south, and the Yanqing District in the northwest will further increase the impact on the surrounding UER. Under the EPS, there will be a slight change in the UER-CDB, but it will generally remain the same as in 2020. It again indicates that the constraints of the EPS can positively maintain the UER level.

4. Discussion

4.1. Implications of Spatiotemporal Variation Characteristics of Land Use

From 2000 to 2020, land use in Beijing was dominated by forest and was mainly located in the western and northern mountainous areas, the area of which increased slightly (Figure 5). This is mainly because the mountainous areas in the west and north of Beijing are unsuitable for large-scale development, limiting uncontrolled land expansion for urban construction. In addition, as ecological problems have become more prominent, the government has designated most mountainous areas as ecological reserves [55]. Implementing various forest protection projects, such as the “Three Norths” protection forest [56], has led to a steady increase in forest cover over the past 20 years. Therefore, in the NDS of 2030, the forest area will continue to grow, but more so in the EPS (Figure 5f).

In the last 20 years, the changes in the building area and cropland have been the most obvious. The encroachment of new buildings on a large area of cropland was the main feature of land use change in Beijing, and it mainly occurred in the plains. This feature was also confirmed by Liu et al. [57] and Wang et al. [55] in their studies on land use change in Beijing. As early as 2014, Zhao et al. [58], by studying the spatial and temporal characteristics of cropland change in China over the past 30 years, found that the resultant decrease in cropland in China’s peri-urban areas was mainly associated with an increase in building. In addition, Chu et al. [59], by predicting the spatial distribution of land use in the Beijing–Tianjin–Hebei region under the natural development state in 2030, also found that urban building expansion had the most severe impact on cropland. Beijing is a vital core city in the Beijing–Tianjin–Hebei region. Its topography and spatial land use distribution are similar to the Beijing–Tianjin–Hebei region. Therefore, the future encroachment of new buildings on cropland under the NDS will be also even more severe in Beijing (Figure 5d). We, therefore, set up the EPS with consideration not only for the protection of ecological space but also to strictly limit the loss of basic agricultural land within the cropland protection red line. We found that under the EPS, the rate of the future expansion of the building area will decrease, but its process will continue, and the impact on cropland remains (Figure 5e). However, compared with the NDS, the building area will be limited to 3670 km², which is more in line with all the requirements of Beijing’s urban plans. Chen et al. [29] also predicted the spatial distribution of land use in Beijing in 2030 under three scenarios—natural development, cropland protection, and ecological protection—and concluded that future land use planning in Beijing should take ecological protection and cropland protection into account. He also believed that Beijing’s future land use planning should take ecological protection and cropland protection into account, which can ensure regional ecological security and improve the quality of cropland and protect food security at the same time. This is consistent with the viewpoint of this paper.

4.2. Impact of Urban Land Use on the UER

In this paper, we constructed the UER based on the landscape pattern index and ecosystem service functions, and we found that the spatial distribution characteristics of UER (Figure 6) were consistent with the spatial distribution pattern of land use (Figure 5), which indicates that UER is strongly influenced by land use. Urbanization is the primary driver of land-use change in urban areas [60]. The critical process of rapid urbanization is the conversion of natural vegetation, such as cropland and forests, into buildings, which inevitably leads to changes in regional landscape patterns [61] and ecosystem service functions [62]. Li et al. [63] analyzed the impact of urbanization on landscape patterns in Beijing and found that with the building area expansion, landscapes with ecological functions, such as cropland and forest, became more fragmented and irregular. Xie et al. [64] also demonstrated the loss of ecosystem services due to the replacement of forest, grass, and cropland with the building due to urban expansion in Beijing, a trend that will continue to increase in the future. These studies also further confirmed a strong relationship between the UER and the spatial distribution and change in land uses.

Our results showed that High and Higher resiliencies were mainly in the western and northern mountainous areas, Moderate and Lower resiliencies were mainly in the border between mountainous and plain areas, while Low resilience was in the urban building areas. Since forest patches in mountainous areas are well preserved [63], ecosystems are more stable. The “source” resilience is much stronger than “sink” resilience, hence the highest UER. The periphery of the central and southeastern plains is dominated by the cross-distribution of building areas and cropland [65], and the ecosystem is more disturbed, resulting in lower UER. Once the area is developed on a large scale, the impact of urbanization on UER will be powerful. The urban center is mainly covered by buildings, with less ecological space, and the “source” resilience is lower than the “sink” resilience, making the ecosystem vulnerable to human activities and resulting in the lowest UER. It shows that Beijing’s UER is spatially heterogeneous under the impact of land use. Shi et al. [66] assessed the ecological resilience of the Beijing–Tianjin–Hebei region and found that it was high in the northwestern part of the region and relatively low in the southeastern part. The finding of this study agrees with our results. The study also demonstrated that ecological resilience is negatively correlated with land expansion in urban building areas. In other regions, a study by Huang et al. [67] also pointed out that the ecological network resilience of the Shandong Peninsula urban agglomeration showed a negative correlation with the level of urbanization. Duo et al. [14] took Nanchang as an example and found that the cold spots of UER were mainly located in urban areas and the hot spots were mainly distributed in areas with low land use intensity. The above findings reconfirm that the spatial distribution pattern of UER in this paper is strongly impacted by the urban land use changes.

In 2030, compared with the NDS, the EPS is more able to maintain the stability of the UER as it protects ecological space such as cropland and forest and limits the further expansion of the building area (Figure 6d,e). Xia et al. [33] simulated the spatial distribution pattern of UER in Hangzhou under different future development scenarios and found that the development scenario with protection conditions has a higher UER than the natural scenario, which is consistent with the findings of this paper. However, what cannot be ignored is that, in the future, rapid global urban expansion will still lead to significant losses of cropland and forested land, with various negative impacts [64,68]. In order to reduce ecological risks, it is necessary to accurately identify the location of cropland and forest vulnerable to urban expansion [69]. In the future, emphasis should be placed on strengthening the protection of mountains, water, lakes, fields, forests, and grasses in urban areas, especially in the plains, by strictly guarding the cropland protection red line, promoting ecological restoration projects and continuously improving the quality of ecological space.

4.3. Differences of UER-CDB in Different Development Directions

This paper analyzed the UER-CDB in different development directions in Beijing using the transect gradient method, and showed that the UER-CDB varies significantly in different development directions and that the “source” and “sink” resiliencies do not strictly follow the urban gradient pattern in some directions. Qiu et al. [70] analyzed the landscape pattern changes in eight directions in Beijing using the transect gradient method, and found that each landscape pattern index did not change in different directions according to the urban gradient pattern, but fluctuated, and pointed out that this phenomenon was related to the transect step size, and the shorter the transect step size the greater the fluctuation. After several experiments, we found that the fluctuations in the “source” and “sink” ductility curves stabilize in most directions for a transect step size of 5 km. However, there were still multiple intersections of the two curves in the northeast, south, southwest, and northwest directions (Figure 7b,e,f,h). Our analysis revealed that this phenomenon is caused by the actual development of the city, and there are three main reasons for this. First, such as the Capital International Airport in the northeast and Daxing International Airport in the south, these areas are located at the edge of the central city, and the construction and operation of the airport have increased the resilience of the “sink” factor in this area, while the resilience of the “source” factor in the surrounding areas remains high. Secondly, the Miyun District in the northeast, Fangshan District in the southwest, Changping District in the northwest, and Yanqing District in the northwest are suburban areas, but their development intensity is greater than that of the surrounding areas, so the resilience of the “sink” factor is higher, but its impact on the resilience of the “source” factor is much less extensive than that of the urban center. Chen et al. [71] explored the relationship between ecological quality and road density in Fuzhou City and found that the critical distance between the two is much larger in the urban center than in the suburban counties (ecological quality and road density can correspond to “source” and “sink” in this paper), which is consistent with the results of this paper. However, in the process of urban expansion, the Changping and Fangshan Districts became gradually connected to the urban center, making the UER-CDB increase and spatially exceed these two districts. The development area of the Miyun District and Yanqing District has gradually expanded, making the “sink” resilience gradually increase and the “source” resilience gradually decrease. Third, in the case of the Wenyu River Park, located 20 km to the northeast, the urban expansion has led to an increase in UER-CDB, and the former UER-CDB has been planned and constructed as an inner-city ecological park, resulting in a high-level of the “source” resilience in the area.

In terms of temporal changes, we found that the UER-CDB in Beijing gradually moved from the central urban areas of the Chaoyang District, Haidian District, and Fengtai District to the new urban development areas of the Changping District, Shunyi District, Tongzhou District, Daxing District, and Fangshan District with the urban expansion over the past 20 years. The increase was especially larger in the northeast, southwest, and northwest (Figure 7b,f,h, Table 4). Hu et al. [72] analyzed the land use changes in Beijing from 2005 to 2015 and found that urban expansion mainly occurred in the northeast, southwest, and northwest of Beijing, which is consistent with the results of our study. In the northeast, the urban building area has gradually expanded to the far suburbs due to the economic pull of the Capital International Airport in the surrounding area. As an important node in the northeastern development belt of Beijing and one of the critical new cities, the Shunyi District has developed more rapidly in the last decade. The building has gradually extended to both sides of the Chaobai River through the Shunyi District [73], which has also impacted the ecological pattern of the surrounding area [74]. The increase in the UER-CDB in the southwest is mainly impacted by the Beijing–Hong Kong–Macao Expressway and the Beijing–Kunming Expressway, which are located on the “Beijing-Baoding-Shijiazhuang” axis and are essential areas for promoting the collaborative development of Beijing, Tianjin, and Hebei [59]. The increase in the UER-CDB in the northwest is mainly impacted by the development of the Changping District. The Changping District, with its high topography in the northwest and low topography in the southeast, has a similar topographic structure to Beijing and is also a key development area in Beijing, with development and construction in its southeast reaching the foot of the mountains and the ecological barrier of the mountains limiting the building area’s further expansion.

Furthermore, our findings also found that the UER-CDB in the eastern and southeastern directions changed significantly over the 20 years, gradually approaching Hebei Province or being adjacent to it (Figure 7c,d, Table 4). The Beijing Urban Master Plan (2004–2020) explicitly focused on the development of the Tongzhou District in the east and Yizhuang New Town in the southeast, and the development of Tongzhou into a sub-center of Beijing, undertaking the decentralization of population and functions from the central urban area and the construction of a new industrial agglomeration. The Beijing Urban Master Plan (2016–2035) further suggests that the sub-city center will be a world-class modern urban area and a model of harmonious coexistence between humans and nature. These policy plans will inevitably promote the urbanization of the area. Liu et al. [75] projected the land use development of the Tongzhou District in 2035 and found that the area and distribution of ecological land use under the sustainable urban development scenario was optimal and highly consistent with the planned land use. This scenario is similar to the EPS in our study. It again confirms that an EPS based on ecological protection policies is more in line with Beijing’s future urban development trend and is more conducive to achieving sustainable urban development.

4.4. Advantages and Limitations

The urban ecosystem is a complex network system composed of diverse natural and human-made elements [28]. Traditional UER assessment methods tend to focus on the own characteristics of ecosystems; For example, Duo et al. [14] established a quantitative measurement and assessment framework for UER based on indices such as habitat quality, stability of ecosystem landscape structure, and ecological vitality coefficients, which were integrated from three dimensions of ecosystem resistance, resilience, and vitality, respectively. Shi et al. [66] directly used net primary productivity (NPP) as the base data for assessing UER. These studies ignored the interaction of various factors in urban ecosystems. However, in this paper, the protection and disturbance factors in the urban ecosystem were considered, the factors in the system were divided into “sources” and “sinks” based on the “source-sink” theory, and finally the UER was assessed based on the interaction between the two. In addition, the transect gradient analysis method could capture the trends of protection and disturbance factors in urban ecosystems in different directions more precisely. By this method, we identified the specific locations of UER-CDB and revealed the spatial heterogeneity and formation mechanisms of urban development and ecological patterns in specific directions. Therefore, compared with the traditional methods, the method in this paper has a higher accuracy and stronger relevance.

However, some things could still be improved. Firstly, there is the question of the input data and parameters of the PLUS model. The existing GlobeLand30 data were able to meet the needs of macro-scale studies with our manual correction, though. However, the accuracy of land use data is the main factor affecting the accuracy of the PLUS model prediction. These data still have a gap with the actual land use at the fine scale, and it is difficult to meet the fine-scale analysis and application. In the future, when available, we will examine these data by using high-precision land cover data, such as the basic geographic national monitoring data and the Third National Land Survey data, to improve the accuracy of model simulations. Additionally, when using the PLUS model to carry out simulations, many parameters are set subjectively, which can also produce errors in the simulation results [29]. For example, we changed and calculated the land use transfer probability matrix predicted by the Markov model and finally obtained the land use demand under the EPS. In fact, many factors affect land use changes in the development of cities, so there is still a certain gap between the simulation results and the actual land use. In the future, we will try to use other models to predict land use demand under different future development scenarios, such as system dynamics (SD) [38], multi-objective planning (MOP) [76], and others, and introduce multiple driver data to improve simulation accuracy.

Secondly, for the UER assessment method, we only based the land use data, although this can reflect the impact of land use and its changes on UER. However, assessing UER through land use change alone may ignore other ecological information. Therefore, in future studies, we consider taking more ecological factors into account, especially multiple ecosystem services functions such as water production services, carbon sequestration services, and food services, in order to better reflect the complexity and diversity of urban ecosystems. Meanwhile, we obtained a time-slice result, so we will try to further consider the time scale and dynamic changes to better capture the changing characteristics of urban ecosystems.

Finally, it has some limitations in identifying UER-CDB. We only identified and analyzed the change characteristics of UER-CDB in eight urban development directions, but this inevitably ignored some details in other directions, making the study results not comprehensive enough. In addition, we used a rectangular quadrat of 5 km × 8 km for the transect gradient analysis, and although the study conditions were satisfied after the experiment, this would inevitably lose some details, so we consider identifying and analyzing UER-CDB at a more comprehensive and finer spatial scale in future studies.

5. Conclusions

This study took Beijing as an example. Firstly, we used the Markov matrix and PLUS model to simulate the spatial distribution of land use in 2030 under the NDS and EPS. Then, based on the “source-sink” landscape theory, the “source” resilience values and “sink” resilience values were calculated by combining the landscape pattern index (basic resilience) and ecosystem service function (process resilience), and then the UER was calculated comprehensively. Finally, the UER-CDB was identified in eight directions using the transect gradient method: east, northeast, north, northwest, west, southwest, south, and southeast. The following conclusions were drawn:

• We demonstrate that the effective assessment of UER and quantitative identification of UER-CDB can be achieved by using the “source-sink” theory and the transect gradient method and successfully applying it in practice to historical periods and different development scenarios in the future. This study is an important complement to the methodological study on the assessment of UER and the identification of its critical thresholds.
• Over the past 20 years, land use in Beijing was dominated by forest, accounting for more than 40% of the total area, and this proportion has been showing an increasing trend. The encroachment of new buildings on cropland has been the main feature of land use change in Beijing, mainly in plain areas. In the next ten years, compared with the NDS, ecological spaces such as cropland, forest, grass, and water will be strictly protected under the EPS, and land expansion for the building will slow down.
• High and Higher UER areas in Beijing are mainly located in the western and northern mountainous areas, but the area share decreased from 63.26% to 48.15% from 2000 to 2020. In contrast, Low UER areas are mainly located in the building areas of the city, with an increased share of 11.94%. In the future, the EPS restrains the expansion of the building area, resulting in UER in Beijing in 2030 remaining stable compared with that in 2020.
• From 2000 to 2020, the changes in UER-CDB in Beijing in different development directions had obvious differences. Among them, the increase in the northeast, southwest, and northwest was more than 25 km. The changes were also obvious in the east and southeast due to the impact of policy planning. In the future, the UER-CDB will further increase under the NDS, while it remains basically the same as in 2020 under the EPS.
• Compared with the NDS, the EPS based on ecological protection policies is more in line with Beijing’s future urban development plans. It has proven to be both ecologically safe for the region and to improve the quality of cropland, which is more conducive to sustainable urban development.