Multiscale Analysis of the Effects of Landscape Pattern on the Trade-Offs and Synergies of Ecosystem Services in Southern Zhejiang Province, China

Abstract

Identifying the trade-offs and synergies (TOSs) of ecosystem services (ESs) and their responses to landscape patterns at various scales, especially in mountainous areas, could benefit the strategies of ES management and landscape optimization. In this study, the southern Zhejiang Province, a hilly region in eastern China, was chosen as the study area. Five ESs, including food production (FP), carbon sequestration (CS), flood mitigation (FM), water conservation (WC), and soil retention (SR) in 2020 were quantified. The TOSs of these ESs were identified at four spatial scales (i.e., grid, watershed, town, and county scales) through Pearson correlation analysis and the spatial overlay method. The effects of landscape patterns on the TOSs of ESs were analyzed by applying a logistic regression model. Results showed that FP and other ESs were trade-offs, while the other ES pairs were synergies. Spatial overlay results showed that weak synergies increased significantly, while strong synergies decreased significantly with the increase of the scale. The direction of the influence of landscape pattern on TOSs did not change, but the magnitudes of the impacts were scale-dependent. Landscape composition (i.e., cropland%, forest%, construction land%) had more significant effects on the trade-offs of ESs than spatial configuration (i.e., LSI, PD, COHE, and SHDI). The magnitudes of impact of landscape composition were strengthened at larger scales, while the effects of landscape configuration on the TOSs of ESs became complex as the scale changed. The results of this study could contribute to understanding how landscape patterns affect TOSs across scales, which will promote the hierarchical governance of ESs in mountainous areas.

1. Introduction

Ecosystem services (ESs) refer to the multiple direct or indirect benefits that people obtain from ecosystems [1,2]. Land use/land cover (LULC), topographic conditions, and climate change play an important role in the formation and allocation of ESs [3]. Inappropriate changes in LULC may lead to the extinction of local species and the reduction of natural habitats and ecosystem functions, thereby affecting the provision of ESs, such as carbon sequestration [4], water conservation [5], and food production [6]. As ES provision is based on the interaction of multiple related ecological processes, complex relationships usually exist among ESs, i.e., trade-offs and synergies (TOSs) [7]. Moreover, landscape patterns have been greatly transformed and fragmented by LULC changes and human activities [8]. Previous studies have shown that landscape pattern can affect ESs and TOSs, and a higher fragmentation or lower connectivity could lead to fragile ecosystems and could impair the formation of ESs [9,10]. It is urgent and critical for decision-makers to have a deep understanding of how landscape patterns affect the TOSs of ESs in order to achieve a win-win outcome through management strategies. Scale effects cannot be ignored when analyzing the relationship between landscape patterns and ecosystems because landscape patterns are scale-dependent and TOSs are spatially heterogeneous [11,12]. Therefore, a multiscale analysis of the effects of landscape patterns on TOSs is essential for ecological management and decision-making [13].

Landscape patterns (e.g., components and configurations) have been identified as one of the most substantial drivers of TOSs [1,14]. Several studies have evaluated the interactions between landscape patterns and ESs [15,16]. Lamy et al. (2016) revealed that landscape composition impacted the spatial pattern of ESs due to the increased competition of ESs [17]. Yang et al. (2018) supported this conclusion, since they found that excessive land reclamation enhanced food production but weakened other services, such as habitat quality and carbon storage [18]. Additionally, studies that focused on landscape configuration agreed that configuration altered ESs as well [14,19]. The effect of landscape patterns on TOSs among ESs has recently been studied, with findings showing that landscape patterns can either enhance or reduce trade-offs/synergies [10,20]. However, spatial heterogeneity has not been sufficiently considered to study how landscape patterns affect TOSs, especially in multiple scales.

Spatial scales are considered an essential factor to understand and practice sustainability in an operational way [21]. Studies focused on associations between ESs at specific scales, such as county scale [22], watershed scale [23], grid scale [24], and plot scale [25], finding that TOSs between the same pair of services vary across different scales. For example, Hou et al. (2017) found that there were trade-offs between water yield, habitat quality, and evapotranspiration at the pixel scale, but they vanished at the town scale [26]. Other studies indicated that ES synergies were enhanced as the scale increased [27,28]. Such a difference mainly stems from the unconformity of ecological conditions and processes at different scales or the spatial mismatch between the supply and demand of ESs [29]. In general, local residents pay more attention to provision services they can directly enjoy, such as food supply, while decision-makers attach more importance to regulating services (e.g., water conservation, flood mitigation, etc.) that are related to long-term well-being [30,31]. In this case, the preferences of different stakeholders should be comprehensively considered in decision-making, and the emphasis and priority of services should be weighed to maximize the benefits of overall ESs. Therefore, scientists should recognize that TOSs can change over space, and it is even more remarkable to consider scale effect when analyzing the influence of landscape patterns on TOSs among ESs.

TOSs among ESs have been evaluated in many studies through statistical analysis, such as the Spearman correlation analysis [32,33,34], Pearson correlation analysis [35,36], and other statistical approaches [37]. Correlation analysis is a widely used method for identifying TOSs, but it can only reflect the linear relationship of each pairwise ES [38]. To measure the relationships of different ESs, researchers used spatial overlay analysis by setting some explicit manners, and they can visually show the spatial patterns of TOSs [39]. For example, Zhang et al. (2022) evaluated the TOSs of forest ecosystem services by using spatial overlay analysis, which helps to weigh management decisions corresponding to specific spatial locations, so as to implement forest resource management more effectively [29]. Wang et al. (2022) demonstrated that TOSs showed obvious spatial differentiation, which depended on the impact of urban development and terrain conditions [40]. Therefore, the spatial overlay method could analyze the relationship of multiple ESs and more effectively demonstrate the spatial differentiation of TOSs, which could provide the basis for decision-makers to weigh the pros and cons of ES management.

Mountainous areas provide various ESs, with rich geographical features and natural resources. Compared to plain ecosystems, mountain ecosystems are more fragile and sensitive and are much more susceptible to environmental changes. Southern Zhejiang Province, which includes 26 mountainous counties (MC), has experienced a fast development since 2001, resulting in dramatic LULC changes [41]. However, there are still ecological and environmental problems, such as soil erosion and forest ecosystem degradation, due to human activities [42]. Therefore, it is essential to further explore the TOSs among ESs in the MC region without reducing ecological and social benefits. Existing studies have focused on quantifying ESs, exploring TOSs, and investigating their driving mechanism [37,43,44]. To our knowledge, to date, few studies have explored the impacts of landscape patterns on ESs or TOSs. To make up for the above limitations, this study aimed to analyzed the relationship between landscape pattern and TOSs at multiple scales and to explain the causes of scale effect. In this study, three research objectives were explored: (1) to examine spatial patterns of multiple ESs; (2) to explore the TOSs at multiple scales; and (3) to analyze the effects of landscape pattern on TOSs at different scales.

2. Study Area and Materials

2.1. Study Area

The MC region lies between 27°06′ N–30°2′ N and 118°1′ E–121°56′ E in the southern Zhejiang Province. It covers an area of 45,707 km², accounting for 43.5% of Zhejiang Province (Figure 1). Zhejiang has promoted the development of less-developed areas as the strategic focus of the province’s modernization. In 2001, the government, for the first time, specified that the less-developed areas included all 16 counties in Quzhou and Lishui City, as well as Taishun, Wencheng, Yongjia, Cangnan, Pan’an, Wuyi, Sanmen, Xianju, Tiantai, and Chun’an counties. In 2005, the objective increased to Pingyang County, namely the 26 mountainous counties.

The study area has a typical subtropical monsoon climate, with an annual average temperature range of 15 to 18 °C and an annual average precipitation of 1100 to 2000 mm. The elevation decreases from southwest to northeast, ranging from 1911 m to 0 m above sea level. Zhejiang’s forest coverage is about 61%. In contrast, the forest coverage in the study area exceeds 76%. The total GDP of the MC region was CNY 591.6 billion in 2020, accounting for 10.6% of the total GDP of the province. The permanent resident population at the end of 2020 was 10.17 million people, accounting for 17% of the province’s total population. The per capita GDP is about CNY 62,100, less than 68.8% of the province’s per capita GDP.

2.2. Data Sources

Data mainly included: (1) The classification of terrestrial ecosystems in Zhejiang Province was based on the 30 m × 30 m land-cover dataset of the Chinese Academy of Sciences (https://www.resdc.cn/) accessed on 20 April 2022. (2) Digital elevation model (DEM) data were derived from the Geospatial Data Cloud, Chinese Academy of Sciences (https://www.gscloud.cn/) accessed on 22 April 2022. (3) Net Primary Production (NPP) was extracted from MODIS/Terra Net Primary Production Gap-filled Yearly L4 Global 500 m SIN Grid (MOD17A3HGF), which was supplied by NASA-USGA (https://lpdaac.usgs.gov/) accessed on 17 May 2022. (4) The normalized difference vegetation index (NDVI) data were extracted from the MODIS/Terra Vegetation Indices 16-Day L3 Global 250 m SIN Grid (MOD13Q1). We calculated annual NDVI data using the maximum value composite method, which was used to generate the vegetation cover factor. (5) Hydrology-related spatial data were produced by the Department of Water Resources of Zhejiang Province. The runoff coefficient for different ecosystems was adopted from the previous study [45]. (6) Meteorological data, including precipitation and evaporation with a spatial resolution of 1 km × 1 km (rainstorm), were obtained from the Meteorological Bureau of Zhejiang Province. (7) Statistical data, such as grain output, were obtained from the Zhejiang Statistical Yearbook.

3. Methods

3.1. Quantifying Multiple ESs

In this study, we identified five ES types based on the following three criteria: (1) The selected ESs should be closely related to human well-being in the MC region. The study area belonged to a water conservation zone in the main function regionalization of China [46], and it was located in a mountainous area, with high carbon sequestration and soil retention. Lakes and reservoirs are densely distributed and have high flood mitigation. We also considered food production services that directly benefit local residents. (2) The data needed to calculate the selected ESs should be available. (3) The selected ESs should be significantly affected by human activities and socio-economic development. Based on this, five key ESs (food production, FP; water conservation, WC; soil retention; carbon sequestration, CS; flood mitigation, FM) were proposed to be selected for use in the MC region.

3.1.1. FP

Food production is a key ES for food security and sustainable development. The food supply in the study area mainly comes from cultivated land, involving grain, beans, oil, and vegetables. We recognized cultivated land as the main land-use type providing food production, and there existed a linear relationship between food production and NPP [47]. We extracted the NPP of each grid by the average value in August, which reflected the best growth status. Then, the food crop production of the whole county was allocated to the farmland grids according to the NPP as follows [48]:

$$F_i=NP{P}_{ij}\times \frac{G_j}{NP{P}_{sum,j}}$$

where F_i is the food production in pixel i (t), NPP_ij is the farmland NPP in pixel i in county j (t), NPP_sum,j is the total value of farmland NPP in county j (t), and G_ij is the grain yield in county j (t).

3.1.2. WC

Water conservation refers to water held by forest, shrub, grassland, and wetland ecosystems from precipitation. We counted the quantity of annual water conservation by using water balance model, and runoff coefficient values were obtained from the previous study [49].

$$Q_{\mathit{wc}}={\displaystyle \sum }_{i=1}^n{A}_i{\times (P}_i-R_i-{\mathit{ET}}_i{)\times 10}^{-3}$$

where Q_wc is the quantity of water (m³·a⁻¹) conserved by an ecosystem, A_i is the area of the ecosystem i (m²), P_i is the precipitation of the ecosystem i (mm·a⁻¹), R_i is the runoff of the ecosystem i (mm·a⁻¹), which is related to precipitation using the runoff coefficient, ET_i is the evapotranspiration of the ecosystem i (mm·a⁻¹), and n is the number of ecosystem types.

3.1.3. SR

Soil retention is a function of the ecosystem that holds soil in place to keep it contained. Due to the large amount of hills and mountains, abundant rainfall, and anthropogenic activity in the MC region, soil erosion has become one of the major problems, resulting in the loss of limited arable land. Therefore, soil retention is an important ecosystem service that should be valued in the study area. We used the Revised Universal Soil Loss Equation (RUSLE), which can be expressed as follows [50]:

$$Q_{\mathit{sr}}={\displaystyle \sum }_{i=1}^{n}R\times K\times L\times S\times \left(1-C\right){\times A}_i$$

where Q_sr represents the amount of soil (t·a⁻¹) held by the ecosystem. R is the rainfall erosivity factor (MJ·mm ·ha⁻¹·h⁻¹·a⁻¹) based on monthly rainfall, K is the soil erodibility factor (t·ha·h·ha⁻¹·MJ⁻¹·mm⁻¹) derived from Zhang et al. [51], L is the slope length factor, S is the steepness factor calculated by ARCGIS 10.2, C is the vegetation cover factor [52], A_i is the area of the ecosystem i (m²), and n is the number of ecosystem types.

3.1.4. CS

Carbon sequestration refers to the long-term capture and storage of atmospheric carbon dioxide in terrestrial ecosystems. NPP is one of the factors that can reflect the total amount of CO₂ that is sequestrated [53]. In this study, we calculated the amount of sequestrated carbon as follows:

$$Q_{{\mathit{CO}}_2}{=M}_{{\mathit{CO}}_2}{/M}_C\times \alpha \times \mathit{NPP}$$

where $Q_{{\mathit{CO}}_2}$ is the quantity of the sequestrated carbon (t·CO₂·a⁻¹) by certain ecosystem, $M_{{\mathit{CO}}_2}/M_C$ is the ratio of CO₂ and C, $\alpha$ is the conversion coefficients of net ecosystem production (NEP) and NPP [54], and NPP is the annual net primary productivity (t·C·a⁻¹).

3.1.5. FM

Flood mitigation refers to the storage capacity of ecosystems, which can store excess water flows during flooding. Ecosystems such as forests, lakes, marshes, and reservoirs can store water and affect the timing and magnitude of water flows and runoff, intercepting heavy rainfall and absorbing water through roots and storage capacity [55]. Regardless of no marsh land cover in this study area, we measured flood mitigation based on the runoff retention by natural vegetation, lakes, and reservoirs [56].

Q_fm = Q_fc + Q_lc + Q_rc

where Q_fm is the annual quantity of flood mitigation provided by the ecosystem (m³·a⁻¹), and Q_fc, Q_lc, and Q_rc are the annual quantity of flood mitigation from natural vegetation, lakes, and reservoirs, respectively (m³·a⁻¹). The raster calculation results of Q_fc, Q_lc, and Q_rc were applied to the raster calculator in ARCGIS 10.2 for raster superposition to obtain Q_fm.

For natural vegetation, the quantity of flood mitigation from natural vegetation was measured as the difference between runoff without and under vegetation cover. We calculated the quantity of flood mitigation as follows:

$$Q_{\mathit{fc}}={\displaystyle {\displaystyle \sum }_{i=1}^n}\left(P_i-R_{\mathit{fi}}\right){\times A}_i{\times 10}^3$$

$$R_{\mathit{fi}}=P_i{\times \beta }_i$$

where P_i is the precipitation of ecosystem i (mm·a⁻¹), R_i is the runoff of ecosystem i (mm·a⁻¹), β_i is the average surface storm runoff coefficient of ecosystem i [57], A_i is the area of ecosystem i (km²), and n is the amount of ecosystem types.

For lakes, we calculated the quantity of flood mitigation from lakes based on the relationship between storage capacity and lake area. For lakes in Zhejiang, the relationship was constructed as follows:

$$Q_{\mathit{lc}}={\displaystyle {\displaystyle \sum }_{i=0}^n}{e}^{4.924}\times A_i^{1.128}\times 3.19$$

where A_i is the area of lake i (km²), and n is the total number of lakes.

For reservoirs, the flood control storage capacity was used as the quantity of flood mitigation from reservoirs.

3.2. Selection and Calculation of Landscape Pattern Metrics

According to previous studies, changes in landscape patterns affect various ecological processes and functions, ultimately leading to changes in ecosystem services [58]. Based on existing studies and combined with the actual situation of the study area [8,9,58,59], seven landscape pattern metrics were selected for analysis, i.e., the proportion of cropland (cropland%), proportion of forest (forest%), proportion of construction land (construction land%), patch density (PD), cohesion (COHE), landscape shape index (LSI), and Shannon’s diversity index (SHDI), respectively. Among them, cropland%, forest%, and construction land% were landscape composition indicators, which represented the dominant type of regional land use. PD, COHE, LSI, and SHDI were landscape configuration indicators. Descriptions of all indicators are provided in

Table 1. Metrics used for quantifying landscape patterns.

Metrics	Description (Unit)	Equation	Description of the Parameter in the Calculation Formula
Cropland%	Proportion of cropland (%)	$\mathrm{Cropland}\%=\frac{Are{a}_{\mathit{Cropland}}}{Are{a}_{total}}$	Areacropland represents the area of cropland; Areatotal represents the total area.
Forest%	Proportion of forest (%)	$\mathrm{Forest}\%=\frac{Are{a}_{\mathit{Forest}}}{Are{a}_{total}}$	Areacropland represents the area of forest; Areatotal represents the total area.
Construction land%	Proportion of construction land (%)	$\mathrm{Construction}\mathrm{land}\%=\frac{Are{a}_{\mathit{Construction}\mathit{land}}}{Are{a}_{total}}$	Areacropland represents the area of construction land; Areatotal represents the total area.
PD	Landscape fragmentation (n/km2)	$\mathrm{PD}=\frac{N}{A}$	N represents the number of landscape patches; A is the total landscape area.
COHE	The connectivity of patches (%)	$\mathrm{COHE}=\left[1-\frac{{{\displaystyle \sum }}_{j=1}^m{P}_{ij}}{{{\displaystyle \sum }}_{j=1}^m{P}_{ij}\sqrt{a_{ij}}}\right]{\left[1-\frac{1}{\sqrt{Z}}\right]}^{-1}\times 100$	Pij is the perimeter of the patch; aij is the area of the patch; Z is the number of cellular.
LSI	The complexity of patch shape (unitless)	$\mathrm{LSI}=\frac{0.25P_{ij}{{\displaystyle \sum }}_{k=1}^m{e}_{ik}^{*}}{\sqrt{a_{ij}}}$	Pij is the perimeter of the patch; e*ik is the edge in the landscape between class i and k; aij is the area of patch.
SHDI	Landscape diversity and the extent to which the landscape is dominated by a few landscape types (unitless)	$\mathrm{SHDI}={\displaystyle {\displaystyle \sum }_{i=1}^m}{P}_i{\times \mathrm{lnP}}_i$	Pi is the proportion of landscape occupied by class i.

.

3.3. Analyzing Trade-Offs and Synergies of ESs at Different Scales

To identify the quantitative correlations between different ESs, Pearson correlation analysis was used to reveal the interactive relationships between ESs. A correlation analysis was conducted at the grid, watershed, town, and county scales, and the mean value of each ES in each unit was determined.

To further reveal the spatial interactive relationships between ESs, the spatial overlay analysis was introduced to reveal the TOSs among the multiple services. First, the value of each ES category was divided into three levels using the natural break classification method, namely low (L = 1), medium (M = 2), and high (H = 3). Then, the raster data of the five ecosystem services were overlaid as follows:

CODE = FM × 10,000 + SR×1000 + WC × 100 + CS × 10 + FP

where FM, SR, WC, CS, and FP represent flood mitigation, soil retention, water conversation, carbon sequestration, and food production, respectively. CODE is a five-digit code, representing the supply capacity of the corresponding ES type.

Table 2. Classification of the spatial TOSs of ESs.

Relationship	Classification	Spatial Combination	Samples
Trade-offs	Strong trade-offs	1H and 4L; 1H, 1M, and 3L; 1H, 2M, and 2L; 1H, 3M, and1L.	11,311; 11,113; 11,321; 12,113; 12,312; 12,321; 22,312; 32,212.
	Weak trade-offs	2H and 3L; 2H, 1M, and 2L; 2H, 2M, and 1L; 3H and 2L; 3H, 1M, and 1L; 4H and 1L.	11,133; 11,313; 23,113; 31,123; 23,213; 23,123; 33,113; 23,313; 33,213; 33,133; 33,313; 33,133.
Synergies	Weak synergies	1M and 4L; 2M and 3L; 3M and 2L; 4M and 1L; 5L.	11,211; 12,111; 12,112; 12,211; 12,212; 22,112; 22,212; 22,122; 11,111.
	Strong synergies	5H; 4H and 1M; 3H and 2M, 2H and 3M; 1H and 4M; 5M.	33,333; 33,233; 23,333; 33,223; 32,332; 22,332; 32,322; 22,322; 22,232; 22,222.

lists the classification criteria for the types of TOSs.

3.4. Logistic Regression Model

An ordered logistic regression model was used to examine the TOSs of ESs and landscape patterns. A logistic regression model is a classical machine learning method. It builds a linear regression based on the Sigmoid function, which is usually used to address the problem of the multi-variable quantitative analysis of binary classification (dependent variable y = 1, 0). The logistic regression model has great advantages over support vector machines and neural network models in terms of training and recognition time. In this study, various landscape metrics (landscape composition indicators and landscape configuration indicators) were used as the independent variables, and the types of TOSs among ESs were applied as the dependent variables. The specific formula was as follows:

$$\begin{gathered} \lg P=Z=C+B_1{X}_1+B_2{X}_2+\dots +B_n{X}_n \\ P=\frac{e^Z}{1+e^Z} \end{gathered}$$

where P is the probability of the types of TOSs, X is the landscape metrics, and B is the logistic regression coefficient.

4. Results

4.1. Spatial Patterns of Multiple ESs

The spatial distribution of five ESs exhibited clear heterogeneities (Figure 2). As shown in Figure 2a, food production was located in the plains in the northwest of the MC region, such as Longyou, Qujiang, and Jiangshan, as well as Pingyang and Cangnan in the eastern coastal plains. Carbon sequestration, flood mitigation, water conservation, and soil retention displayed similar spatial patterns. In detail, the values of carbon sequestration service were higher in the western parts of the study area that were covered by forest, while the values were lower in the central parts of Quzhou and the coastal areas in Cangnan and Pingyang (Figure 2b). As shown in Figure 2d, a higher water conservation value was observed in the region’s northwestern and southern regions, including Chun’an, Kaihua, Changshan, Suichang, Yunhe, Longquan, Qingyuan, Jingning, Taishun, Wencheng, and Cangnan. The high-value regions of soil retention were located in the mountainous areas that were covered with dense forest (Figure 2e).

4.2. Correlation Relationships between ESs

As shown in Figure 3, significant correlations were found among 10 pairs of ESs (Pearson correlation coefficient p < 0.05; −1 < r < 1) at each scale. FP, which was negatively correlated with other ESs (p < 0.05), had a high negative correlation with FM at the grid scale (p < 0.05, r = −0.55). WC, FM, SR, and CS were positively correlated with each other, and FM was positively correlated with CS (p < 0.05, r = 0.92) and WC (p < 0.05, r = 0.87) to a great extent. When expanded to the watershed level, the correlation coefficient between FP and FM increased (p < 0.05, r = −0.49), and the positive correlations between FM and CS (p < 0.05, r = 0.93), as well as WC (p < 0.05, r = 0.88), were enhanced. At the town level, the negative correlation between FP and other ESs was stronger, while the positive correlation between WC, FM, SR, and CS showed different trends. For example, the coefficient value between FM and CS was essentially the same, the coefficient value between CS and SR decreased, and the coefficient value between WC and SR increased. The positive correlations between WC, FM, SR, and CS were all weakened when the scale was increased to the county level.

4.3. TOSs of ESs at Different Scales

Table 3. The percentage of TOSs at the grid, watershed, town, and county scales.

TOSs	Strong Trade-Offs	Weak Trade-Offs	Weak Synergies	Strong Synergies
Grid scale	13.40%	46.03%	20.50%	20.07%
Watershed scale	14.69%	46.05%	19.80%	19.45%
Town scale	13.82%	50.61%	19.94%	15.63%
County scale	17.32%	39.09%	34.35%	9.24%

shows that weak trade-off area had the largest proportion through all scales. At the grid scale, 20.07% of ESs had strong synergies with each other and were scattered in forest ecosystems, as well as in mid-altitude hills. A total of 13.40% of ESs had strong trade-off and were distributed in agricultural and urban ecosystems (Figure 4a, Table 3). There was a total of 19.45% ESs that had strong synergies in the watershed level and were mostly located in mountainous areas of mid-low altitude in the western area and eastern coast (Figure 4b, Table 3). Table 3 and Figure 4c show that the weak trade-off area was mainly distributed in the southwest and northwest of the study region. When expanded to the county level, 34.35% of ESs had low synergies and were mainly distributed in Chun’an, Changshan, Wuyi, Songyang, Liandu, Jinyun, Pan’an, Tiantai, and Sanmen (Figure 4d, Table 3). Meanwhile, 9.24% of ESs had strong synergies and were displayed in the southern part, such as Wencheng, Taishun, and Cangnan.

In order to show the variation of trade-offs and synergies at different spatial scales, we further used a Sankey diagram to clearly represent the increase or decrease of trade-offs and synergies at the grid, watershed, town, and county scales (Figure 5). In general, the TOSs were dominated by weak trade-offs, which first increased from the grid scale to the town scale and then decreased from the town scale to the county scale (Figure 5, Table 3). The increase in weak trade-offs was mainly attributed to the inflow of strong synergies at the town scale. The weak synergies of ESs experienced a decrease from the grid scale to the watershed scale, followed by a significant increase from the town scale to the county scale. Specifically, approximately one-third of the weak synergies of ESs at the county scale came from weak trade-offs at the town scale. Strong synergies showed a decrease from the grid scale to the county scale.

4.4. The Impacts of Landscape Patterns on TOSs at Different Scales

Table 4. Ordered logistic regression results ^a at the town, watershed, and grid scales.

Driving Factor	Strong Trade-Offs		Weak Trade-Offs		Weak Synergies		Strong Synergies
	β	Exp(β)	β	Exp(β)	β	Exp(β)	β	Exp(β)
	Town scale
COHE	–	–	0.523	1.687	−0.762	0.467	–	–
LSI	–	–	–	–	–	–	0.819	2.269
PD	–	–	−0.831	0.435	–	–	–	–
SHDI	–	–	–	–	–	–	–	–
Cropland%	0.900	2.459	−1.349	0.259	–	–	–	–
Forest%	−1.899	0.150	1.024	2.783	–	–	–	–
Construction land%	–	–	–	–	0.492	1.635	−0.681	0.455
	Watershed scale
COHE	1.421	4.141	–	–	−0.609	0.544	–	–
LSI	–	–	–	–	–	–	0.842	2.322
PD	0.607	1.835	–	–	−0.442	0.643	–	–
SHDI	0.954	2.597	−1.443	0.236	1.008	2.741	−0.459	0.632
Cropland%	0.523	1.687	−1.789	0.167	−0.244	0.784	0.281	1.324
Forest%	−0.704	0.495	0.745	2.107	–	–	0.150	1.162
Construction land%	–	–	–	–	0.235	1.265	−0.536	0.585
	Grid scale
COHE	–	–	0.642	1.900	–	–	0.718	2.049
LSI	–	–	−0.694	0.500	0.556	1.744	0.619	1.858
PD	0.221	1.248	–	–	–	–	–	–
SHDI	−0.370	0.691	–	–	0.681	1.975	0.941	2.564
Cropland%	0.521	1.683	−1.833	0.160	−0.458	0.632	0.569	1.766
Forest%	−0.522	0.593	0.556	1.744	−0.530	0.589	1.542	4.672
Construction land%	−0.233	0.792	−0.275	0.759	0.242	1.274	−0.833	0.435

Exp(β) value is used as an effect size indicator to measure the influence of an independent variable on a dependent variable. “–” indicates that driving factors that do not contribute significantly to explaining the landscape matrix were excluded from the final regression equation. ^a All variables (driving forces) are significant at p < 0.01.

shows the regression coefficients of the impacts of landscape patterns at different scales, which could indicate the relationship between landscape patterns and the trade-offs/synergies of ESs. Landscape composition (i.e., cropland%, forest%, and construction land%) had a dominant impact on the trade-offs of ESs compared to spatial configuration (i.e., LSI, PD, COHE, and SHDI).

At the town scale, the results showed that cropland% had a positive effect on strong trade-offs (p < 0.01, β = 0.900), and construction land% had a positive effect on weak synergies (p < 0.01, β = 0.492). Forest% was found to have a negative impact on strong trade-offs (p < 0.01, β = −1.899). In terms of the landscape configuration, PD was negatively correlated, and COHE was positively correlated with weak trade-offs (p < 0.01, β_PD =−0.831, β_COHE = 0.523).

At the watershed scale, the positive impact of cropland% and construction land% on strong trade-offs or weak synergies were both weakened (p < 0.01, β_cropland% = 0.523, β_{construction land%} = 0.235). Furthermore, SHDI was found to have a positive effect on strong trade-offs (p < 0.01, β = 0.954). COHE was positively correlated with strong trade-offs (β = 1.421, p < 0.01).

At the 3000 m grid scale, the negative impact of forest% on strong trade-offs was weakened (p < 0.01, β = −0.522), but its positive impact on strong synergies enhanced obviously (p < 0.01, β = 1.542). Additionally, the negative impact of construction land% on strong synergies increased (p < 0.01, β = −0.833). In terms of landscape configuration, COHE was found to have a positive effect on strong synergies (p < 0.01, β = 0.718). The impacts of PD and SHDI on strong trade-offs were both weakened (p < 0.01, β_PD = 0.221, β_SHDI = −0.37).

5. Discussion

Our research provided a multiscale approach to clarify the effect of landscape pattern on TOSs among ESs at different scales. Our findings enrich existing knowledge on the scale-dependent changes of ESs and TOSs, shedding new light on how landscape patterns affected TOSs at various scales, and it provides insights on integrating scale dependence into the sustainable hierarchical governance of ESs [60,61].

5.1. Multiscale Characteristics of ESs Trade-Offs and Synergies

Scale is a key factor affecting the ES characteristics and TOSs because of the unconformity of ecological conditions and processes at different scales or because of the spatial mismatch between the supply and demand of ESs [7,11,62,63]. Combined with the results of this paper, we found that there were distinct differences in the TOSs of ESs at different scales. The synergies between FM-CS and FM-WC increased when the scale expanded from grid to watershed (Figure 3). The spatial overlay results showed that the weak synergies of ESs experienced an increase from town scale to county scale (Figure 5, Table 1). These results were supported by Yang et al. (2021), who found that the ES synergies enhanced as the scale increased, and this may be attributed to the ‘peak cutting and valley filling’ process of map scale synthesis. In addition, the synergies between CS and SR decreased from watershed level to county level (Figure 3), and the strong synergies of ESs showed a continuous decrease from grid-scale to county-scale (Figure 5, Table 1). This was similar to the study of Su et al. (2020), who found that ES synergies between sediment retention and water yield decreased with the increasing scale. The reason may be due to the gradual impact of socio-economic factors such as agricultural production, population, and non-agricultural output on ESs as the scale increases [63]. The changes in the significance of the relationships among ESs illustrated that the potential TOSs between some ESs could be obscured at larger spatial scales [64,65,66,67].

While correlation analysis can calculate the trade-offs or synergies among ESs, it can hide the spatial heterogeneity of ES’s relationships [68,69]. For example, the correlation results showed that 10 pairs of ES relationships between the five ESs types were robust across scales (Figure 3) mainly due to land-use consistency [13,28,66]. However, the spatial overlay results showed that changing the scale can enhance, weaken, or even reverse the interaction between ESs (Figure 5), which suggested that the relationship of ESs changed both in direction and strength across scales [70,71]. For example, approximately one-third of the weak synergies of ESs at the county scale came from weak trade-offs at the town scale. This illustrated that the optimal management of ecosystems can transform trade-offs at the regional scale into opportunities for synergistic enhancement between ESs at the local scale [72].

The trade-off between ecosystem services is spatially heterogeneous and changes with the passage of spatial scales [73]. Therefore, there is a need to link scale effects with specific geographical environments to provide insights for cross-scale ecosystem management and sustainable spatial planning [72,74,75]. According to our methods, weak synergies meant that at least one service was low (L) and lacking ESs at the high (H) level, which underlined the necessity of land management. As areas with weak synergies were mostly located in a region with a relatively developed economy and intensive human activities, food production and regulating services were both at a relatively low level (Figure 2 and Figure 4). Therefore, for these areas, on the one hand, it is necessary to strengthen the management of land-use types, especially at the county scale, set up ecological red lines [76], implement the Grain for Green Project (GFGP) [77], improve the supply capacity of stock, carbon sequestration, water conservation, and flood mitigation services, and alleviate the strong trade-off between services. On the other hand, the government could improve the utilization rate of water and soil resources through scientific irrigation and farming methods, aiming to achieve the high yield, high quality, and high efficiency of crops.

5.2. Multiscale Analysis of the Effects of Landscape Pattern on the Trade-Offs and Synergies

Analyzing the relationship between the TOSs among ESs and the landscape pattern can be a new approach to improving and coordinating ESs [9]. This study analyzed the effect of landscape patterns on TOSs by setting multiple scales. The results showed that the direction of the influence of landscape pattern on TOSs did not change, but the magnitudes of the impacts were scale-dependent.

The findings of this study showed that landscape composition and spatial configuration had significant impacts on the TOSs of ESs, especially landscape composition (Table 3). This might be because landscape composition affected ESs by directly changing the underlying surface, which increases the competition of ESs [17]. In particular, cropland was positively correlated with the strong trade-offs of ESs (Table 3) because the use of cropland under the pressure of human demand pays more attention to food production, a large amount of fertilizers, pesticides, and high-intensity cropland utilization models, leading to the destruction of cropland ecosystems and resulting in the decline of the WC, SR, and CS of cropland. However, in the ecologically dominant functional area, cropland was negatively correlated with the weak trade-offs of ESs, indicating that cropland may be conducive to alleviate the potential conflict between ESs in ecologically sensitive areas (Figure 4, Table 3). It could be explained by the fact that the cropland in these areas is mainly used by plantations, and the planting of a large number of trees increases vegetation coverage to a certain extent, which can improve the regulation and regulation services [78]. Furthermore, the effect of forest area on the trade-offs of ESs is opposite to that of cropland (Table 3), indicating that the impact of forest and cropland on ESs changed in mountainous areas and was at a trade-off condition, which needs addressed to alleviate the conflict between farming and forestry. Fragmented landscapes can hinder inter-regional energy flows and alter nutrient cycling processes, thereby alleviating ecosystem services [79]. PD and SHDI had a positive impact on the strong trade-offs at the watershed scale. This might be because the expansion of construction land and cropland in the areas with strong ES trade-offs led to the increase of PD and SHDI, which could destroy the natural habitats and reduce the ecosystem functions [80,81].

The magnitudes of impact of landscape composition (i.e., cropland%, forest%, and construction land%) were strengthened at larger scales. Specifically, the influence of cropland% (positive) and forest% (negative) on strong trade-offs increased simultaneously from the grid scale to the town scale, and construction land% had a stronger negative impact on high synergies from the watershed scale to the town scale. The main reason was that landscape composition changed when grid cells were aggregated to large scales due to the ‘peak cutting and valley filling’ process, thus determining the ecological effects at each scale [82,83]. According to Figure 4, 13.40% of ESs had strong trade-offs and were distributed in agricultural and urban ecosystems. Specifically, most of the area was located in the middle and lower reaches of the river, which was suitable for providing food production services. However, due to the relatively gentle topography and frequent human activities, the forest stands have poor quality, mainly consisting of secondary forests and artificial plantations of young ages [42]. Therefore, the supply capacities of the CS, SR, FM, and WC were relatively low, resulting in prominent strong trade-offs. In order to alleviate the strong trade-off among these services, policy makers should strengthen the management of secondary forests and plantation forests at the town scale, promote mixed coniferous and broadleaved forests to replace pure forests, and improve the supply capacity of CS, SR, and WC [29]. Nevertheless, at the 3000-m grid scale, construction land% had the strongest impact on high synergies of ESs. This may be because the areas with high synergies were mainly distributed in forest-dominated areas (Figure 1 and Figure 4), and the construction land expansion in these areas directly led to the reduction and fragmentation of natural habitats, resulting in the degeneration of SR and CS [8,84].

The effect of various landscape configuration factors on the different TOSs of ESs was complex as the scale changed. As the scale rose to the town scale, the negative impact of COHE on weak synergies was strengthened. COHE is an indicator of physical connectivity between patch types, and a higher value indicates better connectivity between patches [21]. Therefore, the advantage of COHE can be maximized for landscape configuration planning at the town scale. When the scale increased from the grid scale to the watershed scale, the impact of LSI on the strong trade-offs of ESs was strengthened. The main reason was that landscape pattern indices changed with different scales, and in the process of small-scale transformation to large-scale, the size and current situation of the study unit changed significantly, and the shape of the landscape became complex, thus causing the change of LSI [27,85]. These findings suggest a more precise and detailed landscape configuration design to better assist the synergetic ecosystem services. According to the spatial heterogeneity of ESs, the sustainable planning strategies of landscapes should be adopted in different scales. Emphasis should be placed on improving grain production services, avoiding disorderly expansion and the fragmented development of construction land, and increasing spatial aggregation [16].

5.3. Limitations and Future Work

This research provided a deep understand of the relationships between landscape patterns and TOSs of ESs under various spatial scales in the southern Zhejiang Province of China. These findings may lead to additional knowledge that will provide more operational options for the hierarchical governance of ESs. For landscape compositions, administrative regions (i.e., town level) should be used as ecosystem service management units. For landscape compositions, landscape sustainable planning strategies should be adopted in administrative regions, as well as in watershed and grid scales. However, there are several limitations that need to be investigated with caution. Identifying the mechanisms driving service trade-offs between landscape patterns, natural processes (e.g., precipitation), and socioeconomic development (e.g., population) can help develop scientifically efficient management solutions. Apart from landscape patterns, other driving mechanisms of TOSs among ESs have not been quantitatively explored in this paper. Therefore, in the future, the extent and differences of the influence of landscape pattern, natural environmental factors, and socioeconomic factors on the TOSs should be explored simultaneously to clarify the main driving factors [86]. Second, multiscale spatial analysis is meaningful to identify the effective spatial scales for ES management [66]. Based on the characteristics of ecosystem landscape pattern at different scales and ecosystem trade-offs/synergies, this paper proposed hierarchical ecosystem governance at four scales according to local conditions. However, in terms of sustainable ecosystem supply and regional sustainable development, a more precise and detailed policy design needs to be made by detailed scale analysis. For example, continuous buffer scale (e.g., 5 km, 10 km, 15 km, etc.) analysis may bring a new insight in the research. Third, our study was based on data collected in 2020. The interpretation of the relationship between ESs, TOSs, landscape patterns, and temporal dynamics may be the focus of future research.

6. Conclusions

In order to develop a deeper understanding of the scale-dependent effect of landscape patterns on trade-offs and synergies among ESs, this study investigated spatial patterns, relationships between ES pairs, TOSs among the multiple ESs, and the effect of landscape pattern on the TOSs of ESs at different spatial scales in the MC region. These findings were fundamentally important for guiding local governments to refine and optimize regional TOSs, especially in developing areas that are facing a series of land-use changes because of dramatic development. The results showed that:

(1) There were distinct differences in the TOSs of ESs at different scales. The Pearson correlation results showed that 10 pairs of ES relationships were robust across scales, while the spatial overlay results showed that TOSs changed both in direction and strength across scales. The weak synergies of ESs experienced a significant increase from the town scale to the county scale, while the strong synergies of ESs showed a continuous decrease from the grid scale to the county scale.

(2) The direction of the influence of landscape pattern on TOSs did not change, but the magnitudes of the impacts were scale-dependent. Landscape composition (i.e., cropland%, forest%, and construction land%) had more significant effects on the trade-offs of ESs than spatial configuration (i.e., LSI, PD, COHE, and SHDI). The magnitudes of impact of landscape composition (i.e., cropland%, forest%, and construction land%) were strengthened at larger scales. The effects of various landscape configuration factors on the TOSs of ESs were complex as the scale changed.

(3) A more accurate and detailed landscape pattern design is proposed for sustainable landscape planning strategies across different scales. Specifically, for landscape pattern compositions, administrative regions (i.e., town level) should be used as ecosystem service management units due to the impact of landscape pattern composition on trade-off synergies at larger scales. For landscape pattern configurations, landscape sustainable planning strategies should be adopted in different scales. Emphasis should be placed on improving grain production services, avoiding disorderly expansion and the fragmented development of construction land, and increasing spatial aggregation.