Spatio-Temporal Changes and Trade-Offs/Synergies among Ecosystem Services in Beijing from 2000 to 2020

: Exploring the dynamic changes and trade-offs/synergies among ecosystem services is essential to urban ecological protection and sustainable development. In this study, we quantiﬁed the spatio-temporal changes in nine ecosystem service values in Beijing from 2000 to 2020 based on land-use data and the equivalent factor method. Correlation analysis and geographically weighted regression were combined to explore the trade-offs and synergies between ecosystem services. The results show that (1) the total ecosystem service value of Beijing increased from CNY 15 billion to 52 billion from 2000 to 2020, and the value was mainly contributed to by forest, cropland, and water. The regulating services covered the largest proportion of the total ecosystem service value, followed by the supporting services; (2) the high-ESV area was mainly located in the mountainous area with abundant forest resources, and the low values were mainly concentrated in central urban areas; (3) most of the ecosystem service pairs had synergies, while the trade-offs mainly existed between food supply services and other services. Measures, such as controlling built-up areas, increasing the area of green space and enhancing vegetation protection, as well as implementing high-quality agriculture, should be taken in order to balance the relationship between ESs and improve ecosystem management in Beijing.


Introduction
Ecosystem service (ES) refers to all kinds of benefits directly and indirectly obtained by humans from ecosystems, which are the basic conditions for maintaining human survival and development in the natural environment [1].However, due to rapid economic development, unlimited human exploitation activities have seriously impacted the natural ecosystem, such as soil erosion, biodiversity loss, and habitat degradation [2].According to the Millennium Ecosystem Assessment (MEA) report, more than half of the world' s ESs are at risk of degradation, with human activities being the most important contributing factor [3].The ecosystem services were classified into four categories, including provisioning, regulating, supporting, and cultural services, in the MEA [3].In 1997, the ecosystem service value (ESV) assessment system proposed by Costanza et al. [1] provided an effective way to assess ecosystem services quantitatively.The accounting techniques of ESV mainly include monetary [4] and non-monetary methods [5].The monetary methods assess the economic value of ESV from the consumer side, while the non-monetary methods evaluate ESV from the production or supply perspective [6].Monetization of the ESV is an effective way to measure an ecosystem' s ability to support economic and social development [7], which can be divided into direct and indirect analyses.Direct calculation methods include the market price method, productivity method, travel cost method, replacement cost method, and willingness survey, which are mainly focused on a specific ES [8].However, it is difficult to obtain data and make a horizontal comparison of the direct methods.In contrast, indirect analyses, such as the equivalent value factor method based on land use/land cover data, Forests 2023, 14, 2314 2 of 18 received more attention due to the relatively easy access to data and low cost [9].Based on the evaluation model proposed by Costanza, Xie et al. proposed an equivalent scale of ESV per unit area of different terrestrial ecosystems in China [10,11], which has been widely used because of its simplicity and reliability, especially for the ESV evaluation results from land use changes.In recent years, the research direction has gradually turned to the impact of ESs on human well-being, and exploring the complex interaction relationship behind ecosystem services has become an important research field.
The complex interaction and mutual restriction among ESs are specifically manifested as trade-offs and synergies; that is, an increase in one ES will inhibit or promote another or more ESs [7].In urban areas, the fragmented distribution of vegetation and water leads to an unbalanced supply of ESs, while the high population density and social-economic activities make it difficult to match the demand and supply of ecological services in space [11,12].A correct understanding of the trade-offs and synergies among ESs is the premise of managing these relationships and is also essential to achieving a balance between the supply and demand of ESs.The existing studies have used various methods, such as correlation analysis [13], spatial autocorrelation and principal component analysis (PCA) [14], and the multi-scenario simulation method [15], to analyze the trade-offs and synergies, among which most studies were conducted on the scale of countries [16], watersheds [17,18], regions [19,20], and urban agglomeration [21].There are few studies at the urban scale, and estimations of ESV are mainly made using static assessments, which lack a dynamic nature in time and space.Beijing, the capital city of China, has been the focus of urban ecosystem research.Some studies have evaluated the ESs [22] and the trade-off/synergy across a transitional area in Beijing [23].However, the analyses of the spatial-temporal characteristics of trade-off and synergy among the ESs in Beijing are still insufficient.The main objectives of this study are (1) to explore the spatio-temporal variations in total ESV and the values of different ESs in Beijing, China, and (2) to analyze the relationship among different ESs qualitatively and quantitatively.We used ArcGIS 10.8 to quantify the spatiotemporal evolution of ESV based on the equivalent factor method and land-use data in Beijing from 2000 to 2020 and then combined correlation analysis and the geographically weighted regression (GWR) method to explore the trade-off and synergy among nine ESs.The research results are expected to improve understanding of the interaction between the ESs and provide a reference for ecological protection in Beijing, China.

Study Area
As China' s political and cultural center, Beijing has 16 districts with a total land area of 16,410.54km 2 .The population of Beijing reached 21.843 million in 2022.Beijing is located in the northwestern part of the North China Plain, surrounded by mountains on the western, northern, and northeastern sides (Figure 1).It has a typical temperate monsoon climate, with an average annual temperature of 12.77 • C, an average annual precipitation of 548.86 mm, and an average annual evaporation of 1785 mm.The Yongding River, the Chaobai River, the North Canal, the Daqing River, and the Jihe Canal form the major river system in Beijing.Forest and cropland are the most typical ecosystems in Beijing.The forest coverage rate reached 44.8% in 2022, and the typical vegetation is deciduous broadleaf forest and warm coniferous forest [24].In recent years, the contrast between rapid social-economic development and the decline in ecological and environmental quality has become increasingly severe.

Data Source and Pre-Processing
In this study, the land-use data from 2000, 2010, and 2020 were derived from the China Land Cover Dataset (CLCD) with a spatial resolution of 30 m × 30 m. CLCD was developed by Yang and Huang [25] based on 335,709 Landsat images and contains annual land cover information of China from 1985 to 2020.There are eight land-use types, including cropland, forest, shrub, grassland, water, barren land, built-up area, and wetland, in the study area.The social and economic data, such as national average grain price (CNY/t), planting area (ha), and annual yield of each grain type (t), were derived from the Beijing Statistical Yearbook and the Data Collection of Cost and Income of National Agricultural Products in corresponding years [26][27][28][29][30][31].We used ArGIS 10.8 to extract the land-use data within the administrative boundary of Beijing, created a fishnet with a grid unit of 1 km × 1 km, and obtained 16,958 grid units.The area of various land-use types in each grid was counted, and we calculated the ESV of each grid unit.

Data Source and Pre-Processing
In this study, the land-use data from 2000, 2010, and 2020 were derived from the China Land Cover Dataset (CLCD) with a spatial resolution of 30 m × 30 m. CLCD was developed by Yang and Huang [25] based on 335,709 Landsat images and contains annual land cover information of China from 1985 to 2020.There are eight land-use types, including cropland, forest, shrub, grassland, water, barren land, built-up area, and wetland, in the study area.The social and economic data, such as national average grain price (CNY/t), planting area (ha), and annual yield of each grain type (t), were derived from the Beijing Statistical Yearbook and the Data Collection of Cost and Income of National Agricultural Products in corresponding years [26][27][28][29][30][31].We used ArGIS 10.8 to extract the land-use data within the administrative boundary of Beijing, created a fishnet with a grid unit of 1 km × 1 km, and obtained 16,958 grid units.The area of various land-use types in each grid was counted, and we calculated the ESV of each grid unit.

Estimation of the ESVs
Costanza et al. [1] proposed the method of equivalent value factor per unit ecosystem area at the global scale, which set the unit ecosystem as the standard functional unit at a large scale that provided ecological service products.The formula for calculating the ESV derived from the model proposed by Costanza et al. [1] is as follows: where ESV and ESV represent the total ESV and the value of the fth ES, respectively.A represents the area of land-use type k.VC and VC are the value coefficients of the ES type f and land-use type k, respectively.The equivalent factor method developed by Xie et al. [10] belongs to the unit-value-

Estimation of the ESVs
Costanza et al. [1] proposed the method of equivalent value factor per unit ecosystem area at the global scale, which set the unit ecosystem as the standard functional unit at a large scale that provided ecological service products.The formula for calculating the ESV derived from the model proposed by Costanza et al. [1] is as follows: where ESV and ESV f represent the total ESV and the value of the fth ES, respectively.A k represents the area of land-use type k.VC k and VC fk are the value coefficients of the ES type f and land-use type k, respectively.The equivalent factor method developed by Xie et al. [10] belongs to the unit-valuebased approach developed by Constanza [1].Based on the framework of Costanza et al. [1] and data from MEA [2], Xie et al. [32] grouped the ESs into four types and nine sub-types for China.In the equivalent factor method, the economic value of each ES from an ecosystem is estimated as the product of equivalent coefficients (dimensionless) and the economic value represented by one standard equivalence factor (CNY•ha −1 •a −1 ).The idea is to establish unified evaluation factors to evaluate the value of ecosystem services.Once the value of any one of the ESVs can be calculated, the other types of ESVs can be calculated accordingly.The equivalent coefficient reflects the relative weight of ESV for a certain ecosystem compared to the standard ecosystem (e.g., cropland) [10,32].Then, the total ESV is summed with the value of different ecosystem services.The overall calculation is as follows: where ESV j indicates the ecosystem service value of ES type j; f jk stands for the ESV equivalent coefficients; S k is the area of land-use type k; E a represents the value of the standard equivalent factor, that is, the economic value of the average annual natural grain yield per hectare of cropland (CNY/ha).
To obtain the equivalent coefficients, Xie et al. [32] conducted a questionnaire survey of more than 700 professionals with an ecological background in China in 2003 and 2008.They developed the table named 'Ecosystem service value per unit area for Chinese terrestrial ecosystems', which shows the equivalent coefficients of all ESs [32].According to the previous study in Beijing [33], we selected nine ESs, including food supply (FS), raw material supply (RMS), air-quality regulation (AQR), climate regulation (CR), regulation of water flows (WFR), soil conservation (SC), nutrient cycling (NC), habitat quality (HQ), and landscape aesthetics (LA), from the four ES groups.The equivalent coefficients for the ecosystem services and land-use types in this study are shown in Table 1.However, this table only concerns the equivalent coefficients at the average level in China, and the ESV is closely related to biomass.When the biomass is greater, the ecological service is stronger.Therefore, we modified the table using the biomass-based regional correction coefficients.According to the regional correction coefficient table by Xie [34] (Table 2), the biomass factor of the cropland ecosystem in Beijing is 1.04.Therefore, the equivalent factor of ecosystem services in Beijing is 1.04 times that of the national level in China.The standard equivalent factor serves as the benchmark for other ESs.The value of the standard equivalent factor of the unit ESV is defined as the economic value of the annual natural grain yield of the national average yield of 1 ha.However, because it is difficult to completely eliminate human interference and accurately measure the food value provided by cropland ecosystems only by natural processes, Xie [10] regards the net profit of grain production per unit area of cropland ecosystems as a standard equivalent factor value, which is equal to 1/7 of the market value of the national per-unit grain yield that year [10].The formula is as follows: E a represents the value of the standard equivalent factor, that is, the economic value of the average annual natural grain yield per hectare of cropland (CNY/ha); i refers to the grain type, p i stands for the national average price of grain type i in a certain year (CNY/t); q i is the yield per unit area of grain type i (t/ha); m i represents the planting area (ha) of grain type i; M stands for the total planting area (ha) for all grain types.
In this study, we used the planting area, annual yield, and national average price of the three main crops of rice, wheat, and maize in Beijing to calculate the economic value.The planting area and annual yield of each grain type (i.e., rice, wheat, and maize) were collected from the Beijing Statistical Yearbook in 2000, 2010, and 2020 [26-28].The national average price of each grain type in the study period was obtained from the Data Collection of Cost and Income of National Agricultural Products [29][30][31].Based on the statistical data and Formula (4), we obtained the value of the standard equivalent factor in Beijing (Table 3).The value coefficient of the ES type and the land-use type per unit area in Beijing were obtained based on the equivalent coefficients, regional correction coefficient, and value of the standard equivalent factor (Table 4).Then, the total ESVs in 2000, 2010, and 2020 were estimated based on the value coefficient and area of each land-use type according to Formula (3).

Ecological Services Change Index (ESCI)
The ESCI describes the change in ES to indicate the relative gain or loss of each ESV [35].The calculation formula is as follows: where

• Static Correlation Analysis
Since the distribution of the geospatial data has non-linear and non-normal characteristics, it is preferred to use non-parametric correlation analysis [36].We used the Spearman correlation analysis to quantify the trade-offs and synergies among nine ESVs.Firstly, the min-max normalization was carried out for nine ESVs, and then sample points were collected based on the 1 km × 1 km fishnet.We used the 'corrplot' package in R 4.2 [37] to conduct the Spearman correlation analysis, which can effectively measure the monotonic relationship between data pairs and capture the non-linear correlation without being sensitive to outliers [38].When the correlation coefficient is positive and passes the significance test, it indicates that the ES pair has significant synergy; when the correlation coefficient is negative and passes the significance test, it indicates that they have a significant trade-off.The greater the absolute value of the correlation coefficient, the higher the degree of trade-off or synergy [39].

• GWR Model
The GWR model was modified based on the traditional regression framework and can test spatial instability [40].In order to quantify the spatial differentiation of the relationship between the ESs, we used the "GWModel" package in R 4.2 [41].Since the ESs are spatially non-stationary and affected by spatial heterogeneity, we used ESVs as the independent and dependent variables in this study.The formula is as follows: where x ij and y i refer to the regression independent variable and the dependent variable; of the jth parameter of point i, reflecting the spatial differentiation of the influence of different parameters.The positive and negative signs of the coefficients represent the positive and negative correlations, and the magnitude shows the correlation degree; (u i , v i ) refers to the spatial location of point i; ε i represents the random error; k is the number of independent variables.

Land-Use Change
Table 5 shows that forest covers the largest area in Beijing, followed by cropland, built-up area, grassland, water, and shrub.The cropland area decreased by 8.4% from 2000 to 2020, while the built-up area increased by 7.9%.The land-use map shows the expansion of the built-up area from 2000 to 2020 (Figure 2).The forest and shrub areas had a minor increase of 1.3% and 0.2%, respectively.The grassland increased by 0.1 from 2000 to 2010 and decreased by 1.2% from 2010 to 2020.The water area decreased by 0.3% and then increased by 0.43% in the second period.The barren land and wetland covered a very small proportion of the study area, and their changes can be disregarded.tionship between the ESs, we used the ʺGWModelʺ package in R 4.2 [41].Since the ESs are spatially non-stationary and affected by spatial heterogeneity, we used ESVs as the independent and dependent variables in this study.The formula is as follows: where x and y refer to the regression independent variable and the dependent variable; β (u , v ) is the intercept constant at point i; β (u , v ) refers to the regression coefficient of the jth parameter of point i, reflecting the spatial differentiation of the influence of different parameters.The positive and negative signs of the coefficients represent the positive and negative correlations, and the magnitude shows the correlation degree; (u , v ) refers to the spatial location of point i; ε represents the random error; k is the number of independent variables.

Land-Use Change
Table 5 shows that forest covers the largest area in Beijing, followed by cropland, built-up area, grassland, water, and shrub.The cropland area decreased by 8.4% from 2000 to 2020, while the built-up area increased by 7.9%.The land-use map shows the expansion of the built-up area from 2000 to 2020 (Figure 2).The forest and shrub areas had a minor increase of 1.3% and 0.2%, respectively.The grassland increased by 0.1 from 2000 to 2010 and decreased by 1.2% from 2010 to 2020.The water area decreased by 0.3% and then increased by 0.43% in the second period.The barren land and wetland covered a very small proportion of the study area, and their changes can be disregarded.According to the Sankey diagram (Figure 3), the transformation directions of the land-use types were similar in the two periods.The decrease in cropland resulted from the conversion into the built-up area (126,671.2ha) and forest (23,265.72 ha).The transformation from cropland into built-up area was the land-use conversion with the largest area proportion, which indicates that a large amount of cropland was occupied in the urbanization process.The increase in forest was mainly due to conversion from the cropland and grassland.Grassland was mainly converted from cropland and forest.
According to the Sankey diagram (Figure 3), the transformation directions of the land-use types were similar in the two periods.The decrease in cropland resulted from the conversion into the built-up area (126,671.2ha) and forest (23,265.72 ha).The transformation from cropland into built-up area was the land-use conversion with the largest area proportion, which indicates that a large amount of cropland was occupied in the urbanization process.The increase in forest was mainly due to conversion from the cropland and grassland.Grassland was mainly converted from cropland and forest.

Temporal Change in the ESVs
From 2000 to 2020, the total ESV in Beijing increased from CNY 15.0 billion to 52.0 billion, with an increase rate of 246.48%, which was mainly contributed by forest, crop, and water (Table A1, Figure 4).Among all the land-use types, forest provided the highest ESV (>70%), followed by water and cropland.The total ESV and the ESVs provided by all land-use types increased from 2000 to 2010.During this period, the ESVs of forest, water, and grassland increased by CNY 56.0 billion, 5.6 billion, and 3.3 billion, respectively.However, the total ESV and the ESVs provided by forest, grassland, water, and shrub decreased during the second period.All the ESVs in 2020 were higher than the values in 2000.

Temporal Change in the ESVs
From 2000 to 2020, the total ESV in Beijing increased from CNY 15.0 billion to 52.0 billion, with an increase rate of 246.48%, which was mainly contributed by forest, crop, and water (Table A1, Figure 4).Among all the land-use types, forest provided the highest ESV (>70%), followed by water and cropland.The total ESV and the ESVs provided by all land-use types increased from 2000 to 2010.During this period, the ESVs of forest, water, and grassland increased by CNY 56.0 billion, 5.6 billion, and 3.3 billion, respectively.However, the total ESV and the ESVs provided by forest, grassland, water, and shrub decreased during the second period.All the ESVs in 2020 were higher than the values in 2000.According to the ESV results (Table A2, Figure 4), the regulating services covered the largest proportion of the total ESV, more than 60%.This was followed by the supporting services (about 25%), while the values of provisioning services and cultural services were relatively low.According to the ESV results (Table A2, Figure 4), the regulating services covered the largest proportion of the total ESV, more than 60%.This was followed by the supporting services (about 25%), while the values of provisioning services and cultural services were relatively low.For the individual ESs, the values of CR, WFR, and SC were higher than the other values.The values of all ESs increased from 2000 to 2010; however, only the value of FS increased from 2010 to 2020.Despite the decrease in most ESVs during the second period, the values of all ESs had an overall increase from 2000 to 2020.The spatial distribution result of the total ESV on the 1 km x 1 km grid (Figure 5 and Table 6) shows that the total ESV in Beijing changed significantly from 2000 to 2020.Each year, the low ESV was mainly located within the urban built-up area, while the high ESV was mostly distributed in the forest, grassland, and water areas.As for the spatial changes in the ESV distribution from 2000 to 2020 (Table 3), the area proportion with the low ESV (<CNY 1000 million) decreased from 54.6% in 2000 to 39.2% in 2020.The area with a relatively low ESV (CNY 1000~3000 million) decreased continuously from 44.4% to 12.5%.The medium ESV (CNY 3000~5000 million) area increased from 0.4% to 46.4%.The area with a relatively high ESV (CNY 5000~7000 million) increased by 14.5% from 2000 to 2010 and then decreased by 14.6% during the second period.There was almost no area with a high ESV (>CNY 7000 million) in 2000.The proportion increased to 35% in 2010 and then decreased to 1.5% in 2020.
According to the ESV results (Table A2, Figure 4), the regulating services covere largest proportion of the total ESV, more than 60%.This was followed by the suppo services (about 25%), while the values of provisioning services and cultural services relatively low.For the individual ESs, the values of CR, WFR, and SC were higher the other values.The values of all ESs increased from 2000 to 2010; however, only the of FS increased from 2010 to 2020.Despite the decrease in most ESVs during the se period, the values of all ESs had an overall increase from 2000 to 2020.

Spatial Distribution of the total ESV
The spatial distribution result of the total ESV on the 1 km x 1 km grid (Figure Table 6) shows that the total ESV in Beijing changed significantly from 2000 to 2020.year, the low ESV was mainly located within the urban built-up area, while the high was mostly distributed in the forest, grassland, and water areas.As for the spatial cha in the ESV distribution from 2000 to 2020 (Table 3), the area proportion with the low (<CNY 1000 million) decreased from 54.6% in 2000 to 39.2% in 2020.The area with a tively low ESV (CNY 1000~3000 million) decreased continuously from 44.4% to 12.5% medium ESV (CNY 3000~5000 million) area increased from 0.4% to 46.4%.The area a relatively high ESV (CNY 5000~7000 million) increased by 14.5% from 2000 to 201 then decreased by 14.6% during the second period.There was almost no area with a ESV (>CNY 7000 million) in 2000.The proportion increased to 35% in 2010 and the creased to 1.5% in 2020.

Distribution of the ESCI
The ESCI distribution in three periods (2000-2010, 2010-2020, 2000-2020) shows the spatial characteristics of gain (ESCI > 0) and loss (ESCI < 0) of the total ESV (Figure 6).During the first period, the ESCI values in 99.2% of the study area were greater than 0, indicating that the total ESV increased in almost all areas.However, the ESCI values in 98.2% of Beijing were lower than 0 from 2010 to 2020, which indicates the general loss of the total ESV.The area where ESV significantly decreased was the built-up area with high expansion intensity.The area with the gain of ESV from 2000 to 2020 covered 95.4% of the study area, and the loss of ESV was mainly distributed in the central urban area and southern part of Beijing.
The ESCI distribution in three periods (2000-2010, 2010-2020, 2000-2020) shows the spatial characteristics of gain (ESCI > 0) and loss (ESCI < 0) of the total ESV (Figure 6).During the first period, the ESCI values in 99.2% of the study area were greater than 0, indicating that the total ESV increased in almost all areas.However, the ESCI values in 98.2% of Beijing were lower than 0 from 2010 to 2020, which indicates the general loss of the total ESV.The area where ESV significantly decreased was the built-up area with high expansion intensity.The area with the gain of ESV from 2000 to 2020 covered 95.4% of the study area, and the loss of ESV was mainly distributed in the central urban area and southern part of Beijing.

Spearman Correlation of the ES Pairs
The Spearman correlation analysis was conducted to quantify the degree and direction of interaction among nine ESs from a static point of view.A total of 36 ES pairs passed the significance test (p < 0.01); that is, all the ESs were significantly correlated.In 2000 and 2010, 28 ES pairs were positively correlated, which indicated that 77.8% of them showed a synergistic relationship (Figure 7).The eight ES pairs related to FS had a trade-off relationship in 2000 and 2010.In 2020, only four pairs were negatively correlated: FS-CR, FS-WFR, FS-HQ, and FS-LA.This indicated that the relationship between FS-RMS, FS-SC, and FS-NC changed from a trade-off to synergy.Among the ES pairs with a trade-off relationship, the correlation coefficients of FS-CR, FS-WFR, FS-HQ, and FS-LA continuously increased, which shows a weaker degree of interaction during the study period.The coefficients of 22 ES pairs with a synergistic relationship increased from 2000 to 2020, indicating a stronger degree of synergy.The coefficients of 10 ES pairs, including RM-AQR, RMS-  The Spearman correlation analysis was conducted to quantify the degree and direction of interaction among nine ESs from a static point of view.A total of 36 ES pairs passed the significance test (p < 0.01); that is, all the ESs were significantly correlated.In 2000 and 2010, 28 ES pairs were positively correlated, which indicated that 77.8% of them showed a synergistic relationship (Figure 7).The eight ES pairs related to FS had a trade-off relationship in 2000 and 2010.In 2020, only four pairs were negatively correlated: FS-CR, FS-WFR, FS-HQ, and FS-LA.This indicated that the relationship between FS-RMS, FS-SC, and FS-NC changed from a trade-off to synergy.Among the ES pairs with a trade-off relationship, the correlation coefficients of FS-CR, FS-WFR, FS-HQ, and FS-LA continuously increased, which shows a weaker degree of interaction during the study period.The coefficients of 22 ES pairs with a synergistic relationship increased from 2000 to 2020, indicating a stronger degree of synergy.The coefficients of 10 ES pairs, including RM-AQR, RMS-NC, AQR-CR, AQR-SC, AQR-NC, AQR-LA, CR-HQ, CR-LA, SC-HQ, and HQ-LA, stayed the same, indicating relatively stable synergies.

Spatial Distribution of Trade-Offs and Synergies
The GWR results show the dynamic trade-offs and synergies among the ESs from 2000 to 2020 (Figure 8).Comparing the proportion of grids with the trade-off and synergistic relationships, the proportion of synergies was higher than trade-offs for most ES pairs except the food supply-related pairs, indicating that most pairs were dominated by the synergistic relationship in space.On the contrary, the ES pairs, including FS-CR, FS-HQ, FS-WFR (in 2000), and FS-LA (in 2000 and 2010), were dominated by the trade-offs, which was consistent with the result of the correlation analysis.The high degree of trade-off grids was mainly distributed in the northern and western mountainous areas of Beijing.As for the ES pairs dominated by trade-offs (e.g., FS-CR, FS-HQ, FS-WFR, and FS-LA), the proportion of trade-offs decreased during the study period, which indicates weaker tradeoffs in 2020.For the ES pairs dominated by synergies, the proportion of synergies increased continuously from 2000 to 2020 (e.g., RMS-AQR, RMS-CR, RMS-WFR, RMS-SC, RMS-NC, AQR-CR, AQR-SC, AQR-NC, CR-WFR, CR-SC, CR-NC, WFR-SC, WFR-NC, WFR-LA, SC-NC, and HQ-LA), indicating that the synergies gradually strengthened.The proportion of synergies for the ES pairs, including RMS-HQ, RMS-LA, AQR-HQ, AQR-LA, CR-HQ, CR-LA, WFR-HQ, SC-HQ, SC-LA, NC-HQ, and NC-LA, increased from 2000 to 2010 and then decreased in the second period although synergies still covered a large proportion.

Spatial Distribution of Trade-Offs and Synergies
The GWR results show the dynamic trade-offs and synergies among the ESs from 2000 to 2020 (Figure 8).Comparing the proportion of grids with the trade-off and synergistic relationships, the proportion of synergies was higher than trade-offs for most ES pairs except the food supply-related pairs, indicating that most pairs were dominated by the synergistic relationship in space.On the contrary, the ES pairs, including FS-CR, FS-HQ, FS-WFR (in 2000), and FS-LA (in 2000 and 2010), were dominated by the trade-offs, which was consistent with the result of the correlation analysis.The high degree of tradeoff grids was mainly distributed in the northern and western mountainous areas of Beijing.A1 and A2).For the change in a single ESV (Table A2), only the value of food supply service underwent a continuous increase from 2000 to 2020.Although the cropland area was smaller in 2020 than in 2000, the value of the food supply service showed an increasing trend between the years.This agrees well with previous findings of increasing food yield in north China in recent decades [42].Changes in the area of cropland are usually the main factor explaining the change in food supply service values [43].However, other important factors also influence food supply service, such as planting density, crop varieties, sowing mechanization, fertilizer application, and so on [44].With limited cropland area, the fertilizer amount in north China increased to support the high grain yields [45], which could explain the increase in the food supply service value despite the decrease in cropland in the study area.The values of other ESs, including raw material supply, air quality regulation, climate regulation, regulation of water flows, soil conservation, nutrient cycling, habitat quality, and landscape aesthetics, experienced a fluctuation; that is, they increased from 2000 to 2010 and then decreased, which was in line with the change in the total ESV.
The change in ESV might be partly due to the change in the value of the standard equivalent factor, which increased in the first period and decreased from 2010 to 2020 (Table 3).From 2000 to 2010, the areas of forest, grassland, and shrub increased by 0.2%, 0.1%, and 0.1%, respectively, while the built-up area expanded by 5.1%.The large increase in the total ESV might result from the increase in the value of the standard equivalent factor rather than the slight increase in natural ecosystem areas.The decrease in ESV in the second period might result from the growth in built-up areas, largely driven by cropland decline.This is demonstrated by many relevant studies.For example, Zhou et al. [46] found continuous expansion of urban built-up areas in cities of the Beijing-Tianjin-Hebei Region from 2000 to 2020, and economic growth promoted urban expansion, especially after the Beijing Olympic Games in 2008.
The spatial distribution of ESV (Figure 5) shows a similar pattern in 2000, 2010, and 2020; that is, high ESVs are mainly concentrated in mountainous and hilly areas in the north and southwest of Beijing, while the low-value areas were located in the central urban area.Xu et al. [47] quantified the habitat quality, water yield, and carbon storage values in Beijing, and the spatial distribution pattern of ESV was similar to that of this study.The high ESV in the mountainous areas might be related to the rich forest resources in these areas.According to the previous study by Chen [48], forest resources show a pattern of more forest in mountainous areas and less forest in plains, as well as more forest outside urban areas and less forest in urban areas.Zhu [49] estimated the forest biomass in Beijing and found that the areas with high biomass were mainly distributed in the north and southwest of Beijing, while the areas with low biomass were mainly distributed in the southeast and central areas of Beijing.Although forests provide more than 70% of the total ecosystem value (Table A1), the forest area remained about 47% of the total area, which did not change much over the study period.The low ESV in the central area might result from the large built-up area and high population density.It has been demonstrated that the specific degree of ESV loss due to the expansion of built-up land is related to the expansion degree of built-up area and the land-use transformation process [50].

Trade-Off and Synergy among the ESs
According to the correlation analysis (Figure 7), there are synergistic relationships among most of the ESs, except the food supply service.Schmalz et al. [51] regard ecosystem services, including water yield, soil retention, water purification, climate regulation, and biodiversity, as water-related ecosystem services (WESs), i.e., products of interactions between water ecosystems and their surrounding terrestrial ecosystems.Food supply is important for food security.It has been found that, when combined with fertilizer application, food production can directly lead to water quality pollution and soil fertility degradation, which can cause permanent damage to WESs [52].The analysis of the relationship between water-related ecosystem services and food supply [53] shows an increasingly prominent trade-off between food supply services and water-related ecosystem services: the food supply caused water pollution and declining soil retention.The balance between food supply and other ESs is an essential prerequisite for achieving both ecological protection and agricultural sustainability.Measures such as the construction of shelterbelts on cropland could be used to reduce the impacts on the other ESs [54].It is also important to increase the level of forest coverage in plains areas in Beijing and use a combination of needle/broadleaf trees and shrubs to improve the ecological benefits.
In all ES pairs with synergies, the correlation coefficient of ES pairs, including the regulation of water flows (WFR), was slightly lower than the coefficients of other ES pairs.The regulation service of water flow is closely related to the vegetation type, climate (e.g., precipitation, temperature), and transpiration [55].The water capture and throughflow control components of flow regulation can be maintained when land cover is modified, provided sufficient vegetation canopy or basal cover is retained.The scarcity of water resources in Beijing could also be a reason for the limitation of water-regulation services.Therefore, increasing precipitation resources, water body area, and vegetation cover can enhance the synergies between regulating water flows and other ESs.
Different ES pairs have different synergies and trade-offs at the grid level (Figure 8).In this study, the type and intensity of relationships among ESs showed spatial heterogeneity (Figure 8).In previous studies, it was found that human activities and ecosystem conditions are the main factors that determine where trade-offs or synergies may occur [56,57].For example, Figure 8 shows a strong trade-off between food supply services and the regulation of water flows, as well as nutrient cycling services in the north and west of Beijing, which is mostly a mountainous area with forests.It has been found that ESs such as food supply, raw material supply, regulation of water flows, and nutrient cycling in suburban areas are influenced by demands for food, recreation, and housing [58].In recent years, many ecological projects, such as the Beijing Plains Afforestation Project [59], the Beijing-Tianjin Sand Source Control Project [60], and the ecological restoration projects on groundwater in the Yongding River Basin [61], have been implemented in Beijing and surrounding cities.However, the trade-offs between ecosystem services were not fully considered in these projects, which caused the incomplete improvement of ESs in some places.In addition to the ecological projects, it is essential to conduct sustainable land-use practices to maximize the benefit of ESs [62].Policymakers should take measures such as controlling the urban built-up area, increasing the green space area and enhancing vegetation protection, implementing high-quality eco-agriculture, and reducing pollutant emissions to balance the relationship between ESs.In areas with ES synergies, policymakers should focus on ecological protection and natural restoration to strengthen the synergy intensity.

Future Outlook
There are still uncertainties concerning the results of this study due to certain limitations.Firstly, the complexity, dynamics, and non-linear characteristics of the ecosystem may introduce defects into the equivalent factor method.Although we used the biomass factor to modify the equivalent coefficient, the local differences should be quantified using more specific parameters.Therefore, a more accurate evaluation model should be developed to improve the comprehensiveness and scientific nature of ESV evaluation in further research.Secondly, this study shows the spatial heterogeneity of trade-offs and synergy.In future research, we will evaluate the ESs through a more comprehensive index system and explore the influencing factors and driving mechanism of the trade-offs and synergies from the perspective of economic, social, and ecological coupling systems.

Conclusions
In this study, we explored the spatial and temporal changes in the ESVs, as well as their trade-offs and synergies, in Beijing from 2000 to 2020.The results show the following: (1) The total ESV in Beijing increased from CNY 15.0 billion to 52.0 billion during the study period, experiencing the process of first rising and then falling.Among all the land-use types, forest provides the highest ESV, followed by water and cropland.The regulating services covered the largest proportion of the total ESV, followed by the supporting services; (2) The spatial distribution of the ESVs in the study area was closely related to land-use types.The highest ESV was distributed in areas with abundant forest resources, and the low ESV was mainly concentrated in urban built-up areas.The area where the total ESV significantly decreased from 2010 to 2020 was the built-up area with high expansion intensity; (3) The static correlation analysis and GWR model indicate that synergy was the dominant relationship between the ESs during the study period, and trade-offs mainly existed between FS and other services.Local food production activities should pay attention to the protection and restoration of the other ESs.The degree and direction of interaction between various ESs changed from 2000 to 2020, and the synergistic degree of most ES pairs strengthened.
These results can help to identify the changing characteristics of various ESs, clarify the interaction process of different ecosystems at the city scale, and provide a reference for targeted policies and measures, such as ecological space regulation, in order to promote ecological protection in Beijing.

Figure 3 .
Figure 3.The Sankey diagram of land-use conversion in Beijing from 2000 to 2020.

Figure 3 .
Figure 3.The Sankey diagram of land-use conversion in Beijing from 2000 to 2020.

Figure 4 .
Figure 4.The percentage of ESVs of each land-use type in Beijing from 2000 to 2020 (FS: food supply, RMS: raw material supply, AQR: air-quality regulation, CR: climate regulation, WFR: regulation of water flows, SC: soil conservation, NC: nutrient cycling, HQ: habitat quality, LA: landscape aesthetics).
For the individual ESs, the values of CR, WFR, and SC were higher than the other values.The values of all ESs increased from 2000 to 2010; however, only the value of FS increased from 2010 to 2020.Despite the decrease in most ESVs during the second period, the values of all ESs had an overall increase from 2000 to 2020.

Figure 4 .
Figure 4.The percentage of ESVs of each land-use type in Beijing from 2000 to 2020 (FS: food supply, RMS: raw material supply, AQR: air-quality regulation, CR: climate regulation, WFR: regulation of water flows, SC: soil conservation, NC: nutrient cycling, HQ: habitat quality, LA: landscape aesthetics).
Spatio-Temporal Changes in the ESVsThe total ESV increased from CNY 15.0 billion to 82.3 billion from 2000 to 2010 and decreased to CNY 52.0 billion in 2020 (Tables

Table 2 .
[34]regional correction coefficients of the equivalent coefficients based on biomass factor in China[34].

Table 3 .
Statistical data used for calculation of the value of the standard equivalent factor [26-31].
ESCI x represents a single ecosystem service change index; ES CUR x stands for the ESV in the final state; ES HIS x stands for the ESV in the initial state.An ESCI value of 0 indicates there is no change in the ESV; that is, no gain or loss.A negative ESCI value indicates the loss of ESV, whereas a positive value indicates the gain in ESV.

Table 4 .
The value coefficient of the ES type and the land use type per unit area in Beijing from 2000 to 2020 (unit: CNY•ha −1 •a −1 ).

Table 5 .
Area and proportion of land-use types in Beijing from 2000 to 2020.

Table 5 .
Area and proportion of land-use types in Beijing from 2000 to 2020.

Table 6 .
Changes in ESV area at different grades from 2000 to 2020.