Nonlinear Effects of Landscape Patterns on Ecosystem Services at Multiple Scales Based on Gradient Boosting Decision Tree Models

Abstract

Exploring the complex effects of landscape patterns on ecosystem services (ESs) has become increasingly important in offering scientific support for effective spatial planning and ecosystem management. However, there is a particular lack of research on the nonlinear effects of landscape patterns on ESs and scale dependence. Taking Huainan City (in China) as a case study, this study adopted the InVEST model to estimate four key ESs including carbon storage (CS), habitat quality (HQ), nitrogen export (NE), and water yield (WY). Then, we calculated the selected landscape metrics at multiple spatial scales. Furthermore, the gradient boosting decision tree (GBDT) model was developed to investigate the relative importance of landscape metrics in explaining ESs and their nonlinear interrelation. The results indicated that most of the selected landscape metrics were significantly correlated with ESs. The GBDT model, which can explore nonlinear relationships, performed better than the linear regression model in explaining the variations in ESs. The landscape-level metrics of the Shannon’s diversity index (SHDI) and the contagion index (CONTAG) and the class-level metrics of the aggregation index (AI) and edge density (ED) were the most important variables that influenced ESs. The landscape metrics affected ESs within a certain range, and the nonlinear effects varied with scale.

1. Introduction

Ecosystem services (ESs), the direct and indirect benefits that human beings receive from ecosystems, significantly impact human well-being and sustainable development [1,2]. With socioeconomic development and the enhancement of human activities, the regional natural ecosystem is increasingly affected by human disturbance [3,4]. According to the Millennium Ecosystem Assessment (MEA) report released by the United Nations in 2005, more than 60% of ecosystem services have been degraded in the past half century [5]. The intensifying conflict between the degradation of ES provision and the continuous increase in human demand results in great pressure on ecosystems and ecological security.

Climate change, soil properties, topography, and land use change have been widely proven to significantly influence the spatiotemporal changes in ESs in various ways [6]. To provide practical suggestions for effective ecosystem management and spatial planning, more attention has been focused on manageable factors. The dramatic changes in landscape patterns due to intense human activities significantly have impacted the structures and processes of ecosystems by modifying the movements of matter and organisms, thereby affecting the provision of ESs [7,8]. Landscape pattern changes are believed to pose great threats to regional ESs. Therefore, the relationship between landscape patterns and ESs can be widely recognized as an effective tool to address sustainable challenges [9,10].

In addition to the impact of landscape composition on ESs, existing studies have proven that ESs are significantly correlated with landscape configuration [11,12]. Landscape fragmentation and connectivity are significantly related to ES changes by varying the shape and spatial arrangement of landscape patches [7,13]. The complex shape of non-arable land patches may reduce soil erosion [14]. Forest complexity and fragmentation are negatively correlated with water retention [15]. An increase in aggregated forestland results in an increase in water yield [16]. In summary, a basic consensus can be reached from previous studies: landscape patterns can significantly influence the distribution and provision of ESs. Various methods have been applied to investigate the relationships between landscape patterns and ESs based on a linear assumption including correlation analysis [17], multiple linear regression [18,19], stepwise regression [16], redundancy analysis [20], and geographically weighted regression [21]. There are limitations of these statistical methods for detecting the impact of landscape patterns on ESs due to simple linear assumptions. Under the impact of landscape pattern changes, different types of ESs show various change dynamics, resulting in nonlinear ES variations [22,23]. A threshold effect exists in the nonlinear changes in ESs [24]. Small changes in landscape patterns could result in significant changes in ESs when thresholds are reached [25]. The nonlinearity of threshold in the impacts are unclear when using statistical methods with linear assumption [8,24]. However, a knowledge gap still exists regarding the nonlinear impact and threshold effect on the ecological effects of landscape patterns, and this information is of great importance for formulating adaptive management policy and spatial planning [23].

In recent years, algorithms such as decision trees [26], random forests [27], and convolutional neural networks [28] have been used to overcome the shortcomings of statistical methods in exploring nonlinear relationships. Using machine learning algorithms that consider nonlinearity has improved our understanding of the relationships between variables [29]. The gradient boosting decision tree (GBDT) model, which is based on decision trees, has recently been used to explain the nonlinear relationship between variables, and this method shows a higher accuracy than random forests [17]. The main advantage over commonly used linear regression models is that the GBDT model can address multicollinearity issues by involving all important explanatory variables [30]. Furthermore, it can reflect the relative contribution of independent variables. In addition, a partial dependence plot generated by the GBDT model allows for visualization of the nonlinear relationships between explanatory and dependent variables, which can improve our understanding of threshold effects [31,32].

The complexity and uncertainty behind the relationships between ESs and landscape patterns are also reflected in the scale dependency because ESs are generated by various ecological processes and structures with different spatial scales [33,34]. The directions and magnitudes of the landscape pattern impacts on ESs are scale-specific and often vary over distinct spatial scales [35]. To date, the relationships between landscape patterns and ESs have usually been explored at a single spatial scale such as the grid, basin, and administrative region scales; however, this approach neglects the dynamics of landscape pattern-driven changes in ESs across various scales [13,16,36]. The scale effects involved in the relationships between landscape patterns and ESs have not been fully addressed [33,34,37]. ESs usually differ across spatial scales as landscape patterns may have different correlations with ESs at different scales [33,38,39]. How ESs respond to landscape patterns across different spatial scales remain poorly understood. The effective management of ESs and landscape planning requires appropriate scale information related to the impact of landscape patterns on ESs [23,40]. Thus, to facilitate scientific landscape planning, the spatial scale needs to be fully considered when investigating the complex relationships between landscape patterns and ESs.

In this paper, we took the typical resource-based city of Huainan in China as a study area, which has experienced significant landscape pattern changes and ES degradation caused by coal mining and socioeconomic development since the 1980s. To fill the above research gaps, the study attempted to comprehensively explore the nonlinear relationship between landscape patterns and ESs as well as analyze their scale effect. The specific objectives of this study were as follows: (1) reveal the spatial characteristics of landscape patterns and ESs in Huainan; (2) explore the nonlinear response of ESs to landscape patterns and its threshold effect; (3) analyze the scale effect of the relationship between ESs and landscape patterns; (4) identify the implications for effective ecosystem management and spatial planning.

2. Materials and Methods

2.1. Study Area

Huainan City (31°54′N–33°00′N, 116°21′E–117°12′E), located in the northcentral part of Anhui Province, eastern China, covers a total area of 2582 km² (Figure 1). The study area is situated in the subtropical monsoon climate zone, with an average annual temperature of 16.6 °C and an annual precipitation of 893.4 mm. The topography of Huainan can be divided into two zones: the northern plain zone with an average elevation of 20–24 m, and the southern hilly zone with an elevation ranging from 17 to 218 m. As shown in Figure 1, Huainan City subordinates the districts of Tianjiaan, Datong, Bagongshan, Xiejiaji, and Panji and the county of Fengtai. The area proportions of cropland, built-up land, waterbody, forestland, grassland, and unused land accounted for 70.90%, 17.65%, 8.84%, 1.48%, 1.07%, and 0.07%, respectively, of the total area in 2020. The population increased from 2.06 million in 2000 to 2.51 million in 2020. The gross domestic product (GDP) increased from RMB 13.28 billion to 111.42 billion, with an average increase rate of 11.22% [41].

Huainan is one of 14 largest coal bases in China. It is rich in coal resources, with a proven reserve of 18 billion tons, accounting for approximately 32% of the reserves in eastern China. Mining may result in serious ecological degradation such as changes in landscape patterns, water pollution, air pollution, soil degradation, land subsidence, and the destruction of biodiversity. By 2018, the subsidence area driven by underground coal mining in Huainan was 298.6 km², and approximately half of this area was submerged in water [42]. Additionally, the area of built-up land increased from 350.82 km² to 455.73 km² and the urbanization rate increased from 44.78% to 61.08% from 2000 to 2020. Long-term intensive coal mining activities and rapid urbanization and development have dramatically influenced the structure and process of regional ecosystems. The contradiction between ecological security and socioeconomic development has become a main obstacle that hinders the sustainable development of Huainan.

2.2. Data Sources

The data from multiple sources used in this study were as follows (Table 1). (1) Land cover data on Huainan for the year 2020 with a spatial resolution of 30 m, obtained from the Data Center for Resources and Environmental Sciences, Chinese Academy of Sciences (RESDC) “http://www.resdc.cn/Default.aspx (accessed on 12 October 2022)”. The land cover data, with an overall accuracy of 88.7%, were interpreted from Landsat 8 remote sensing images. (2) Digital elevation model (DEM) data with a spatial resolution of 12.5 m were obtained from ALOS “https://alos-pasco.com/ (accessed on 12 October 2022)”. (3) Meteorological data for 2020 were obtained from the Chinese National Meteorological Science Data Service Center “http://data.cma.cn/ (accessed on 12 October 2022)”, and the data included the annual average temperature and precipitation. (4) Soil type data with a spatial resolution of 1 km including topsoil sand, silt, clay, and organic carbon fractions were obtained from the Harmonized World Soil Database (HWSD) “https://www.fao.org/home/en/ (accessed on 12 October 2022)”. (5) Socioeconomic statistical data of Huainan including population and GDP for 2020 were obtained from the Huainan Statistical Yearbook. All spatial data were converted to the spatial coordinate system (WGS_1984, UTM Zone 50N), and all raster data were resampled to a resolution of 30 m.

2.3. Quantification of ESs

The key types of ESs were identified based on four criteria [10,20]: (1) strong relevance to human well-being in Huainan; (2) significant relevance to human activities and landscape pattern change; (3) the concerns of local stakeholders regarding their ecological and environmental interests; (4) availability of the data used for quantifying and mapping ESs. Since human disturbance and intense coal mining activities have caused some ecological issues such as water pollution, vegetation deterioration, and biodiversity loss, the vegetation deterioration would threaten the capacity for carbon sequestration. Large-scale subsidence waterlogging caused by underground mining would influence the capacity of water yield. With consideration of the above criteria and the important ecological issues in this study, four ESs were chosen to represent the key services of Huainan City including carbon storage (CS), habitat quality (HQ), nitrogen export (NE), and water yield (WY). Against the background of continuous rapid socioeconomic development, the selected ESs are significantly related to human well-being and play key roles in maintaining ecological security in Huainan. The details of the calculations are listed in Supplementary Section S1. The main paramaters for calculating ESs are presented in Tables S1–S3. To analyze the scale effect of the relationship between landscape patterns and ESs, different scales (grid sizes of 1 km, 2 km, and 3 km) were adopted to map the ESs and the results were further aggregated by using the Zonal Statistics Table tool in ArcGIS 10.5.

2.4. Calculation of Landscape Metrics

Landscape metrics can reflect landscape pattern information. In this study, the landscape metrics were selected to: (1) include a reasonable mix of metrics to characterize the landscape structure; (2) include recommended metrics based on the previous studies [8]; (3) prioritize easily understood and computable landscape metrics; (4) use the smallest number of independent metrics for reducing information redundancy [18,43]. Therefore, five metrics at the landscape level including patch density (PD), mean shape index (SHAPE_MN), contagion (CONTAG), Shannon’s diversity index (SHDI), and aggregation index (AI) as well as four metrics at the class level including edge density (ED), PD, SHAPE_MN, and AI were chosen to quantify the landscape pattern in Huainan City, as these metrics can reveal the complexity, fragmentation, diversity, and connectivity of the landscape (metrics with “_L” represent landscape-level metrics; metrics with “_C” represent class-level metrics). Grids with sizes of 1 km, 2 km, and 3 km were first created by using the Create Fishnet tool in ArcGIS 10.5. The land cover raster data were split into small pieces based on grids of 1 km, 2 km, and 3 km. The selected landscape metrics were further calculated at different scales using Fragstats 4.2.

2.5. GBDT Model

The GBDT model is a recently developed machine learning algorithm that combines the decision tree flowchart approach with the boosting ensemble technique to create a strong classifier from a number of weak classifiers [32]. Compared with traditional linear regression, which has commonly been used in exploring correlations, the GBDT model has many advantages. The GBDT model can simultaneously reveal the relative importance of landscape patterns in affecting ESs and avoid linear assumptions, thus exploring the potential threshold effects of landscape patterns on ESs. In addition, the GBDT model can effectively assess the nonlinear relationships between landscape patterns and ESs. A variation inflation factor (VIF) test should be conducted to ensure that the independent variables are free from multicollinearity if the linear regression model is used. Therefore, the GBDT model can address the multicollinearity issue that usually exists in ES and landscape pattern analysis [44]. In this study, the GBDT model was selected to explore the nonlinear relationship between ESs and landscape patterns.

The GBDT algorithm aims to minimize the loss function $L\left(y,F\left(x\right)\right)={\left(y-F\left(x\right)\right)}^2$ through the function $F_M\left(x\right)$ based on the training sets. In the GBDT model, $F_M\left(x\right)$ can be expressed as follows:

(1)$F_M\left(x\right)={\sum }_{m=1}^M{f}_m\left(x\right)={\sum }_{m=1}^M{\theta }_{m}t\left(x;{\mu }_m\right)$

where M represents the number of decision trees and iterations; $x$ is a set of explanatory variables (i.e., PD, ED, SHAPE_MN, AI, CONTAG, and SHDI); ${\mu }_m$ represents the mean value of split positions and terminal nodes in each tree t(x;μ_m); ${\theta }_m$ represents the weight for evaluating each decision tree and can be estimated by minimizing the loss function.

The gradient boosting approach was developed and applied to estimate the parameters in the above function by reducing the loss function $L$. This method consists of several iterative steps for optimization. The decision tree was initialized as follows:

(2)$f_0\left(x\right)=\beta =argmi{n}_{\gamma }{\sum }_{i=1}^{N}L\left(y_i,\gamma \right),$

The following steps were conducted for iteration round (m = 1, 2, …, M):

Step 1: For each iteration, the negative gradient value of the loss function for data sample $i$ was calculated:

(3)${\gamma }_{im}=-{\left[\partial L\left(y_i,f\left(x_i\right)\right)/\partial f\left(x_i\right)\right]}_{f=f_{m-1}}$

Step 2: The optimal step length θ_m for the search direction of the gradient update was calculated based on the following equation:

(4)${\theta }_m={\mathrm{argmin}}_{\theta }{\sum }_{i=1}^{N}L\left(y_i,F_{m-1}\left(x_i\right)+\theta t\left(x;{\mu }_m\right)\right)$

where ${\theta }_m$ is expected to minimize the loss function value.

Step 3: According to Equations (2)–(4), the regression function was updated:

(5)$F_m\left(x\right)=F_{m-1}\left(x\right)+{\theta }_{m}t\left(x;{\mu }_m\right)$

The “gbm” package in the R platform was applied to develop the GBDT model in this study. A fivefold cross-validation procedure was applied to obtain the optimal parameter settings of the GBDT model. The R² and root mean square error (RMSE) were used to evaluate and compare the performance of the models.

The GBDT model could reveal the relative importance of all explanatory variables. The squared relative importance ($I_{x_i}^2$) of each explanatory variable ($x_i$) was obtained by averaging the square importance value of all of the trees $T$:

(6)$I_{x_i}^2=1/M\times {\sum }_{m=1}^M{I}_{x_i}^2\left(T_m\right)$

The GBDT model can also be used to generate partial dependence plots (PDPs), allowing us to illustrate how an explanatory variable affects the dependent variable. The partial dependence of $F\left(x\right)$ on $xs$ can be expressed as follows [32]:

(7)$F_{xs}\left(xs\right)=1/N\times {\sum }_{i=1}^{N}F\left(xs,x_c^i\right)$

where $xs$ is the explanatory variable that affects the dependent variable, and $x_c^i$ represents other explanatory variables involved in the model. N represents the number of data samples. The “pdp” package in R was applied to generate PDPs, which estimate the average marginal effects of landscape patterns on ESs.

3. Results

3.1. Spatial Patterns of ESs

As indicated by Figure 2, there was significant spatial heterogeneity among the selected ESs in 2020. Higher CS values were predominantly distributed in the hilly region with abundant vegetation coverage, while lower CS values were found in the plain region with a higher coverage of built-up land and water bodies. The hilly areas located in the southern areas had higher HQ values, while the plains areas had lower values. Hilly areas have abundant forests and rich biodiversity, and less land is disturbed by human activities, implying that the ecological structure is less threatened. Lower NE values reflect higher water purification capacity. Higher NE values were observed in cropland because intensive agricultural activity resulted in an increase in pollution exports. The hilly areas with higher vegetation coverage exhibited lower NE values due to the strong water purification capacity of vegetation. Higher WY values were observed in the built-up area because of its lower infiltration capacity compared to the hilly area.

3.2. Spatial Characteristics of Landscape Patterns

Table 2 shows the differences in the spatial metric values among land cover types. The built-up land in Huainan had the highest PD value, which suggests that the fragmentation degree of built-up land was the highest among the different land cover types. The ED and SHAPE_MN values of cropland were the highest, which suggests that the complexity degree of cropland was the highest. The AI value of cropland was higher than that of the other land cover types, indicating the best aggregation. The complexity of the landscape patch shape can be measured by SHAPE_MN. As shown in Table 2, the SHAPE_MN value of cropland was the highest in 2020, implying that the path shape of cropland was the most complex compared with the other land cover types.

The spatial distributions of the landscape metric values at the landscape and class levels for different spatial scales (1 km, 2 km, and 3 km) are presented in Figures S1 and S2, respectively. At the landscape level, the high-value areas of PD were concentrated in the southern part of Huainan, which was occupied by various land cover types. The high-value areas for SHAPE_MN were found in the city center of Huainan, which could be due to the intensive human activity. At the class level, the high-value areas for AI were mainly concentrated in the city center of Huainan. The high-value areas for ED and PD were mainly distributed in the southern part of Huainan. The high-value areas for SHAPE_MN were mainly located around the city center. The statistics of landscape metrics at different scales are presented in Table S4.

3.3. Comparison of the GBDT and Linear Regression Models

To evaluate the performance of the GBDT model, both the linear regression model and the GBDT model were developed based on the landscape metrics and ESs. These models were compared in terms of the R² and RMSE values. A higher adjusted R² and a lower RMSE mean that more of the variance of a dependent variable can be explained by independent variables. As shown in Table 3, the GBDT models showed higher adjusted R² values and lower RMSE values for each spatial scale when compared with the linear regression models, which indicated that the GBDT models performed better than the linear regression models in terms of exploring the relationships between the landscape patterns and ESs.

3.4. Relative Importance of Explanatory Variables

The relative importance values for each explanatory variable were retrieved to represent the impact of the independent variable on the dependent variable. Figure 3 shows the relative importance values of the selected landscape metrics at multiple spatial scales derived by the GBDT model. Although we established three GBDT models for each type of ES involving the same landscape metrics, the landscape metrics made different contributions to each ES in terms of the relative importance values, and the relative importance values varied over different spatial scales. In addition, the relative importance of each landscape metric differed among the four ESs, as indicated by Figure 3. Based on the relative importance values of the landscape metrics, the important landscape metrics were identified.

Compared with the other landscape pattern metrics, SHDI_L was the most valuable variable affecting CS, accounting for 34.19%, 37.44%, and 20.38% at the spatial scales of 1 km, 2 km, and 3 km, respectively. Compared with the other two scales, ED_C exhibited a higher contribution of more than 20% at the 2 km scale. AI_C was also an important variable affecting CS and contributed to 25.90%, 12.51%. and 14.69% of CS at the different scales, respectively. At the 2 km and 3 km scales, SHDI_L (37.49%) and ED_C (22.60%) were the first and second most important variables that impacted HQ, respectively. However, the ranks of these two variables varied over scales, with ED_C becoming the first most important variable at the 2 km and 3 km scales. AI_L contributed more than 6% to HQ at the 2 km and 3 km scales. The greatest contributor to HQ varied from SHDI (1 km) to ED_C (2 km and 3 km). The results indicate that SHDI_L had the greatest impact on NE at the three different scales, accounting for 28.62%, 24.45%, and 31.02%, respectively. The ED_C, SHAPE_L, and SHAPE_C variables contributed more than 10% to NE. AI_C was the most important variable influencing WY, accounting for 38.12%, 36.02%, and 36.91% at the different scales, respectively, followed by ED_C, accounting for 18.24%, 23.34%, and 26.63%, respectively. These results illustrate that the AI and ED of built-up land had a high correlation with WY. CONTAG_L was the third most important variable affecting WY at the 1 km and 2 km scales, accounting for 16.57% and 11.28%, respectively. SHDI_L (6.53%) and SHAPE_L (7.66%) also played important roles in influencing WY at 3 km. The relative importance indicated that the landscape metrics of SHDI, AI_C, ED_C, and CONTAG_L could explain the high variance in the selected ESs. The impacts of PD_C and PD_L were smaller than those of the other variables.

3.5. The Impacts of Landscape Patterns on ESs

The direction and significance of correlations between different ESs and landscape patterns were derived using Spearman correlation analysis and compared at various scales (1 km, 2 km, 3 km). Figure 4 shows that most of the selected landscape metrics were significantly correlated with the four ESs at different scales. Specifically, CONTAG_L was positively correlated with CS, HQ, and NE with relatively higher correlation coefficients, and negatively correlated with WY. However, AI_C had the opposite correlation direction. ED_C had significantly positive correlations with WY. SHDI_L was negatively correlated with CS, with higher correlation coefficients. Furthermore, the results suggested that the direction and magnitude of the correlation between landscape patterns and ESs varied across scales. The strength of the correlations between AI_L and CS decreased as the scale changed from 1 km to 2 km. The correlation between SHDI_L and HQ was also scale dependent. The correlation shifted from negative to nonsignificant and positive along different scales.

The PDPs were further estimated to reveal the nonlinear effects of the three most important landscape metrics on CS, HQ, NE, and WY (Figure 5, Figure 6, Figure 7 and Figure 8). The vertical axis represents the value of each ES, and the horizontal axis represents the value of each explanatory variable. The impacts of landscape patterns on ESs were not linear, and each metric had a significant impact within a certain range. The impacts of the same landscape metric on an ES differed over scale. The results also demonstrated that ESs exhibited heterogeneous relationships with the selected landscape metrics.

Figure 5 indicates that an approximately negative relationship existed between AI_C and CS, reaching a threshold that differed with the spatial scale. SHDI_L exhibited a “U”-shaped correlation with CS. When SHDI_L increased, CS remained unchanged in the range of 0–0.4, decreased in the range of 0.4–0.8, and significantly increased in the range of 0.8–1.0. At the 1 km scale, CONTAG_L was positively correlated with CS, and the correlation was similar to an S-shaped curve, and CONTAG_L had a positive correlation with CS when the value was greater than 40. At the 3 km scale, however, CONTAG_L exhibited significant effects on CS, with rapid increases and decreases in CS when CONTAG_L increased. Both SHAPE_MN_L, and SHAPE_MN_C had an “S”-shaped relationship with CS. At the 2 km scale, for instance, when SHAPE_MN_L ranged from 1.6 to 1.7, a steep positive effect occurred; however, beyond this range, SHAPE_MN_L had minimal effects on CS.

Figure 6 shows the partial dependence of landscape pattern metrics on the HQ. AI_C was negatively correlated with HQ, which was similar to the correlation between AI_C and CS. At the 3 km scale, HQ remained high and steady when AI_C was less than 92. Then, HQ decreased with an increase in AI_C between 92 and 98. HQ remained stable when the HQ exceeded 98. An S-shaped curve was found in the correlation between SHDI_L and HQ. Unlike their correlation at the 3 km scale, SHDI_L exhibited a more complicated relationship with HQ. At the 1 km scale, HQ increased with increasing SHDI_L from 0.0 to 0.3, remained steady when SHDI_L was between 0.3 and 0.7, increased rapidly when SHDI_L was between 0.7 and 0.9, and remained steady beyond this range. Figure 6 indicates that ED_C was negatively correlated with HQ. At the 1 km scale, the increase in ED_C ranging from 1 to 28 was found to have a significantly negative impact on HQ, and the influence was relatively stable within other value ranges. However, the value ranges of ED_C with negative impacts were 4–20 and 10–16 at the 2 km and 3 km scales, respectively.

As shown in Figure 7, SHDI_L exhibited a negative correlation with NE. Specifically, the increase in value ranging from 0.5 to 1.0 at the 1 km scale, the increase in value ranging from 0.6 to 1.1, and the increase in value ranging from 0.7 to 0.9 had significantly negative impacts on NE. The impact on NE was relatively stable beyond this range. Figure 7 suggests that ED_C had a positive correlation with NE. The increases in ED_C ranging from 4 to 30 at the 2 km scale and from 8 to 22 at the 3 km scale were found to have a positive impact on NE. SHAPE_MN_L exhibited a “U”-shaped correlation with NE, which first decreased and then increased. At the 3 km scale, CONTAG_L was found to be positively correlated with NE when the value was in the range of 50–70, and NE remained relatively stable when the value was beyond the range.

Figure 8 displays the relationship between WY and the landscape pattern metrics, which exhibited a significant nonlinear relationship. WY first remained stable and then increased rapidly when AI_C increased. However, the fitted curve varied over scales. In addition, a positive effect existed between ED_C and WY. Specifically, WY increased when ED_C increased within 10–30 (1 km), 5–25 (2 km), and 8–25 (3 km). WY remained constant when ED_C exceeded these ranges. The PDP indicated that CONTAG_L was negatively associated with WY within 30–60 (1 km) and 40–58 (2 km). The impact of CONTAG_L on WY was relatively stable in the other value ranges at the 1 km scale. At the 2 km scale, however, a positive effect of CONTAG_L on WY was found when the value ranged from 70 to 85.

4. Discussion

4.1. The Nonlinear Impacts of Landscape Patterns on ESs

Rapid urbanization, agricultural activity, coal mining, and ecological protection and restoration have been proven to significantly affect landscape patterns [45] and eventually influence ESs [16]. Our results imply that landscape patterns at the landscape and class levels are significantly correlated with ESs, which is in line with the results of previous studies [17,19]. In addition to the effect of landscape composition metrics demonstrated by many studies, relative importance analysis indicated that landscape configuration such as landscape diversity (measured by SHDI), aggregation (measured by AI), connectivity (measured by CONTAG), and complexity (measured by SHAPE_MN and ED) had larger impacts on ESs. In addition to landscape composition, landscape configuration should be fully considered and incorporated when formulating ecosystem management policies and urban spatial planning. The provision of ESs is strongly dependent on the flow of organisms and materials across landscapes, which are affected by landscape patterns [7,19]. The increase in landscape connectivity is beneficial for material transition and ecological processes, further improving ESs generated by ecological processes with actual carriers [46]. The diversity of landscape fragments would increase with the number of classes or the proportional distribution of distinct classes [47]. The mixed configuration of natural and semi-natural landscapes would improve biodiversity and nutrient interception, as indicated by the general increase in habitat quality and a decrease in nitrogen export. The complexity of the patch shape, measured by ED and SHAPE, is sensitive to human disturbance [48]. Human dominated landscapes (i.e., built-up areas) exhibit straight and distinct boundaries, whereas forest land and grassland are usually characterized by irregular shapes and edges [49]. It was further confirmed that different types of landscape configuration accounted for various proportions of the variation in ESs. For instance, the relative contributions of SHAPE_MN and PD for HQ were lower, while the relative contributions of landscape diversity were higher. This result differed from the discovery that HQ was significantly correlated with the landscape complexity and diversity [50]. The magnitudes of the impacts of different variables on ESs could be clearly identified using relative importance analysis.

Furthermore, this study made a unique and significant contribution to investigating the relationships between landscape patterns and ESs through the development of the GBDT model, providing a better understanding of the complex relationships between landscape patterns and ESs. The verification and comparison results confirmed that the GBDT model could be used to explore nonlinear relationships more accurately than other regression models with linear assumptions, which corresponded with previous findings related to GBDT models [31,32,51]. The application of the GBDT model demonstrated that this model not only effectively captured the contributions of different explanatory variables to dependent variables, but also allowed for the visualization of nonlinear effects. Some positive correlations derived from the linear regression models cannot always result in a promoting effect when the variable is manipulated. For instance, landscape diversity measured by the SHDI was positively or negatively correlated with ESs, which was revealed by a linear regression model in previous studies [13,18]. However, the results derived from the GBDT model implied that landscape diversity had a negative effect on CS first and then became a positive effect on CS after reaching a peak value (Figure 5). This is probably because a diverse landscape produces a more complex flow of ecological processes than a landscape with few elements [21]. In addition, our findings were in accordance with previous studies, in which ED was negatively correlated with HQ [16]. The GBDT model further identified the threshold of the response of ED to HQ; beyond this threshold, ED had a minimal effect on HQ. The threshold of ES responses to landscape patterns reflected the phase characteristics of the ecological effects of landscape patterns [23], which need to be fully considered in formulating appropriate urban planning. Based on these findings, a proper threshold for urban planning and design can be identified to mitigate the negative impact of landscape patterns, which was another contribution of this study in assessing nonlinear relationships. Therefore, our findings based on the GBDT model can provide insights into the nonlinear relationships between landscape patterns and ESs, thus supplementing the traditional linear relationship.

4.2. Scale Effect of the Relationships between Landscape Patterns and ESs

Most previous studies have focused on exploring the relationships between landscape patterns and ESs at a single spatial scale [16,52]. Some scholars have argued that the scale dependence effect widely exists in the evaluation of ES and landscape pattern relationships [7,53]. Single-scale analysis may distort the interaction mechanism between landscape patterns and ESs. The results were consistent with previous studies that showed the significant scale effects involved in ESs [52,53,54]. In addition, this study enriches the previous findings that demonstrated the scale dependence of the relationship between landscape patterns and ESs. As indicated by the results, the relative importance of the landscape pattern metrics and their impacts on ESs varied across scales, which demonstrates the sensitivity of landscape patterns and ESs to scale. For example, the correlation between ED and CS changed from significant to nonsignificant at different scales. Although the direction of the response of ESs to some landscape metrics showed good consistency and stability at multiple spatial scales, as shown in the PDPs, the threshold involved in the relationship was different. The comparison of the effects of landscape patterns on ESs among various scales indicated that regulating landscape patterns to promote ESs at fine and coarse scales requires different planning and management policies. If the management policy remains constant across different scales, measures to optimize ESs at one certain scale could lead to a suboptimal result at another scale [7,55]. A mismatch could reduce the efficiency of land planning and management. How to avoid low management efficiency caused by a scale mismatch needs to be fully considered in practice. Our findings imply that understanding the scale effect of the relationships between landscape patterns and ESs will be crucial for identifying suitable scales at which landscape patterns should be regulated to promote and optimize ESs.

4.3. Policy Implications

Effective landscape management is increasingly incorporating urban planning for maintaining or promoting ESs [56]. Given the increasing threat to ESs in Huainan, decision makers aim to maximize ESs by formulating targeted policies and measures. Our findings imply that managing landscape patterns according to relationships at different scales could help optimize ESs, which is similar to the results of previous studies [7,13]. Although relevant information can be obtained based on the results derived from the linear regression model, developing future development policy with support of the GBDT model is potentially more valuable because of the nonlinear effects and thresholds. The exploration of the nonlinear effects of landscape patterns on ESs can prevent the inappropriate interpretation of the impact of landscape patterns on ESs and provide beneficial information for urban planning and management. For instance, the nonlinear relationship between CS and the three most important landscape metrics can be used to determine the suitable landscape metric values for improving CS. In contrast, landscape complexity had a strong negative correlation with HQ within a certain range. Therefore, landscape complexity should be decreased within the threshold range.

The relationships between landscape patterns and ESs strongly depend on the ecological structure and natural environmental conditions. Thus, flexible urban planning and management policies that are adapted to the local conditions are essential for better coordinating the conflict between socioeconomic development and ecological protection. With a large forest coverage area, the southern part of Huainan has higher ecological quality. However, the region suffers from serious challenges from rapid urban expansion. Urban expansion should be strictly controlled by establishing natural preservation zones. In the central part of Huainan, the land was mainly occupied by built-up land with reduced provisioning and regulating services. For highly urbanized regions, some measures such as intensive land use and green infrastructure construction should be implemented to promote ESs. The cropland and coal mining sites are mainly located in the northern part of Huainan. The cropland and ecological land should be reasonably allocated during planning. The “ecological red line”, “permanent basic farmland”, and “urban development boundary” regions should be scientifically identified in planning to realize the coordination among food production, ecological protection, and socioeconomic development. For the coal mining areas in Huainan, the combination of protective exploitation, artificial induction, and natural restoration should be applied as ecological restoration measures in coal mining areas. The subsidence areas caused by underground coal mining need to be reclaimed into ecological parks to further improve the ecological quality of coal mining areas.

4.4. Limitations and Future Work

The GBDT model provides a more reasonable way to investigate the nonlinear correlation between landscape patterns and ESs without assuming a prior function. However, some limitations exist in this study. Some important ESs (e.g., air purification, soil conservation, recreation services) were not considered in this study due to data availability. To provide a comprehensive understanding of ES distribution characteristics and influence mechanisms of landscape patterns, it is necessary to consider more ESs in the future. In addition, the long-term monitoring and analysis of landscape patterns and ESs must be conducted to reveal the temporal dynamics of ESs and their responses to landscape patterns. Moreover, there is a lack of research on the spatial heterogeneity of the relationships between landscape patterns and ESs. The nonlinear effect of landscape patterns on ESs could vary spatially. Therefore, future research should focus on the spatially varying nonlinear effects of landscape patterns on ESs to provide more comprehensive implications for landscape planning and ecosystem management.

5. Conclusions

Understanding the complex relationship between landscape patterns and ESs is vital to effectively formulate and implement spatial planning and ecosystem management. This study applied the GBDT model to analyze the effect of landscape patterns on ecosystem services in Huainan, China. Moreover, this study is only one step toward exploring the nonlinear correlation between landscape patterns and ESs based on a machine-learning model. The GBDT model can overcome the shortcomings of linear regression models and contribute to research on the response of ESs to landscape patterns. We conclude that the performance of the GBDT model is better than that of the traditional linear regression model in explaining the variance in ESs and that applying the GBDT model to assess the nonlinear relationships is more appropriate for interpretation.

The relative importance of each landscape pattern metric was compared to weigh the contribution in explaining ESs. Landscape diversity, complexity, and connectivity were identified as the most important variables that influenced the ESs in the study area. In addition, the nonlinear impacts and potential thresholds of the landscape pattern metrics were visualized through PDPs. SHDI_L exhibited a “U”-shaped correlation with CS. An approximately negative relationship existed between AI_C and CS, HQ, and NE. ED_C had a positive correlation with NE and WY. CONTAG_L had a negative correlation with WY within a certain range. Furthermore, the responses of ESs to landscape patterns varied over scales. Our findings not only demonstrate the usefulness of the GBDT model in exploring nonlinear relationships and threshold impacts but also help to fully understand how landscape patterns influence ESs and provide precise and effective policy implications for decision makers. Understanding the processes that link landscape patterns to ESs will help advance our ability to manage human-dominated landscapes for multiple ESs.

Author Contributions

Conceptualization, C.L. and W.H.; Methodology, C.L.; Software, J.Z.; Validation, C.L. and W.H.; Formal analysis, J.Z.; Investigation, C.L.; Data curation, J.Z.; Writing—original draft preparation, C.L.; Writing—review and editing, W.H.; Supervision, W.H. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Not applicable.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Study area: (a) The location of Anhui Province in the eastern part of China. (b) The location of Huainan City in Anhui Province. (c) The topography and administrative boundary of Huainan City.

Figure 2. Spatial distribution of the four key ecosystem services in Huainan in 2020. (a) The spatial distribution of carbon storage (CS). (b) The spatial distribution of habitat quality (HQ). (c) The spatial distribution of nitrogen export (NE). (d) The spatial distribution of water yield (WY).

Figure 3. Relative importance values of selected landscape metrics at the landscape level and class level for ecosystem services. Metrics with the suffix ‘_L’ represent landscape pattern metrics at the landscape level. Metrics with the suffix ‘_C’ represent landscape pattern metrics at the class level.

Figure 4. Correlation matrix between ESs and landscape pattern metrics at the landscape level and class level at different scales in Huainan.

Figure 5. The partial dependence plots of CS on the three most important landscape metrics at different spatial scales.

Figure 6. The partial dependence plot of HQ on the three most important landscape metrics at different spatial scales.

Figure 7. The partial dependence plots of NE on the three most important landscape metrics at different spatial scales.

Figure 8. The partial dependence plots of WY on the three most important landscape metrics at different spatial scales.

Table 1

Data descriptions and sources used in this study.

Data	Type	Resolution	Data Source
Land cover data	Raster	30 m	Data Center for Resources and Environmental Sciences, Chinese Academy of Sciences (RESDC) “http://www.resdc.cn/Default.aspx (accessed on 12 October 2022)”
Digital elevation model (DEM)	Raster	12.5 m	ALOS “https://alos-pasco.com/ (accessed on 12 October 2022)”
Meteorological data (annual average temperature and precipitation)	Raster	1 km	Chinese National Meteorological Science Data Service Center “http://data.cma.cn/ (accessed on 12 October 2022)”
Soil type data	Raster	1 km	Harmonized World Soil Databased (HWSD) “https://www.fao.org/home/en/ (accessed on 12 October 2022)”
Population and GDP data	Numeric	—	Huainan Statistical Yearbook 2021

Table 2

Landscape metric values at the class level.

Types	PD	ED	SHAPE_MN	AI
Cropland	0.017	19.6782	2.7604	97.8972
Forestland	0.0074	0.6940	1.7807	96.6633
Grassland	0.0089	0.6240	1.7712	95.8261
Waterbody	0.0426	4.0284	1.7978	96.4201
Built-up land	0.6971	15.8244	1.2676	93.3932

Table 3

Comparison of the GBDT model and linear regression model.

Scale	Model	Metric	CS	HQ	NE	WY
1 km	Linear	R2	0.2951	0.5860	0.1716	0.4062
		RMSE	17.1376(Mg/hm2)	0.0706(Unitless)	0.1483(kg/hm2)	87.1089(mm)
	GBDT	R2	0.4986	0.7627	0.3018	0.4859
		RMSE	14.4936(Mg/hm2)	0.0538(Unitless)	0.1356(kg/hm2)	81.9808(mm)
2 km	Linear	R2	0.2117	0.6295	0.3305	0.4462
		RMSE	16.1699(Mg/hm2)	0.0669(Unitless)	0.1030(kg/hm2)	71.8694(mm)
	GBDT	R2	0.4436	0.7635	0.4237	0.6368
		RMSE	13.9704(Mg/hm2)	0.0506(Unitless)	0.0864(kg/hm2)	58.7452(mm)
3 km	Linear	R2	0.2150	0.7357	0.4392	0.4383
		RMSE	13.7860(Mg/hm2)	0.0488(Unitless)	0.0740(kg/hm2)	66.0749(mm)
	GBDT	R2	0.3442	0.8280	0.5285	0.5864
		RMSE	12.5798(Mg/hm2)	0.0437(Unitless)	0.0634(kg/hm2)	55.3898(mm)

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/rs15071919/s1, Figure S1: Landscape metrics at the landscape level at the 1 km scale (left column), 2 km scale (central column), and 3 km scale (right column); Figure S2: Landscape metrics at the class level at the 1 km scale (left column), 2 km scale (central column), and 3 km scale (right column); Table S1: Carbon density per unit area of different land cover types in Huainan City, China; Table S2: Biophysical table in the InVEST NDR model for Huainan City, China; Table S3: Parameters for habitat quality in Huainan City, China; Table S4: The statistics of landscape metrics at different scales [23,57,58,59,60,61,62,63].

References

1. Costanza, R.; d’Arge, R.; de Groot, R.; Farber, S.; Grasso, M.; Hannon, B.; Limburg, K.; Naeem, S.; O’Neill, R.V.; Paruelo, J. et al. The value of the world’s ecosystem services and natural capital. Nature; 1997; 387, pp. 253-260. [DOI: https://dx.doi.org/10.1038/387253a0]

2. Daily, G.C. Nature’s Services: Societal Dependence on Natural Ecosystems; Island Press: Washington, DC, USA, 1997; pp. 454-464.

3. Ding, T.; Chen, J.; Fang, Z.; Chen, J. Assessment of coordinative relationship between comprehensive ecosystem service and urbanization: A case study of Yangtze River Delta urban Agglomerations, China. Ecol. Indic.; 2021; 133, 108454. [DOI: https://dx.doi.org/10.1016/j.ecolind.2021.108454]

4. Vitousek, P.M.; Mooney, H.A.; Lubchenco, J.; Melillo, J.M. Human domination of earth’s ecosystems. Science; 1997; 277, pp. 494-499. [DOI: https://dx.doi.org/10.1126/science.277.5325.494]

5. Millennium Ecosystem Assessment (MEA). Ecosystems and Human Well-Being: Synthesis; Island Press: Washington, DC, USA, 2005.

6. Cord, A.F.; Bartkowski, B.; Beckmann, M.; Dittrich, A.; Hermans-Neumann, K.; Kaim, A.; Lienhoop, N.; Locher-Krause, K.; Priess, J.; Schröter-Schlaack, C. et al. Towards systematic analyses of ecosystem service trade-offs and synergies: Main concepts, methods and the road ahead. Ecosyst. Serv.; 2017; 28, pp. 264-272. [DOI: https://dx.doi.org/10.1016/j.ecoser.2017.07.012]

7. Mitchell, M.G.E.; Bennett, E.M.; Gonzalez, A. Strong and nonlinear effects of fragmentation on ecosystem service provision at multiple scales. Environ. Res. Lett.; 2015; 10, 094014. [DOI: https://dx.doi.org/10.1088/1748-9326/10/9/094014]

8. Lyu, R.; Zhao, W.; Tian, X.; Zhang, J. Non-linearity impacts of landscape pattern on ecosystem services and their trade-offs: A case study in the City Belt along the Yellow River in Ningxia, China. Ecol. Indic.; 2022; 136, 108608. [DOI: https://dx.doi.org/10.1016/j.ecolind.2022.108608]

9. Qiu, J.; Turner, M.G. Importance of landscape heterogeneity in sustaining hydrologic ecosystem services in an agricultural watershed. Ecosphere; 2015; 6, art229. [DOI: https://dx.doi.org/10.1890/ES15-00312.1]

10. Yohannes, H.; Soromessa, T.; Argaw, M.; Dewan, A. Impact of landscape pattern changes on hydrological ecosystem services in the Beressa watershed of the Blue Nile Basin in Ethiopia. Sci. Total Environ.; 2021; 793, 148559. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2021.148559]

11. Eigenbrod, F. 2016. Redefining landscape structure for ecosystem services. Curr. Landsc. Ecol. Rep.; 2016; 1, pp. 80-86. [DOI: https://dx.doi.org/10.1007/s40823-016-0010-0]

12. Shrestha, R.P.; Schmidt-Vogt, D.; Gnanavelrajah, N. Relating plant diversity to biomass and soil erosion in a cultivated landscape of the eastern seaboard region of Thailand. Appl. Geogr.; 2010; 30, pp. 606-617. [DOI: https://dx.doi.org/10.1016/j.apgeog.2010.01.005]

13. Guo, S.; Wu, C.; Wang, Y.; Qiu, G.; Zhu, D.; Niu, Q.; Liu, Q. Threshold effect of ecosystem services in response to climate change, human activity and landscape pattern in the upper and middle Yellow River of China. Ecol. Indic.; 2022; 136, 108603. [DOI: https://dx.doi.org/10.1016/j.ecolind.2022.108603]

14. Jordan, G.; Rompaey, A.V.; Szilassi, P.; Csillag, G.; Mannaerts, C.; Woldai, T. Historical land use changes and their impact on sediment fluxes in the Balaton basin (Hungary). Agric. Ecosyst. Environ.; 2005; 108, pp. 119-133. [DOI: https://dx.doi.org/10.1016/j.agee.2005.01.013]

15. Kuriqi, A.; Pinheiro, A.N.; Sordo-Ward, A.; Garrote, L. Influence of hydrologically based environmental flow methods on flow alteration and energy production in a run-of-river hydropower plant. J. Clean. Prod.; 2019; 232, pp. 1028-1042. [DOI: https://dx.doi.org/10.1016/j.jclepro.2019.05.358]

16. Hao, R.; Yu, D.; Liu, Y.; Liu, Y.; Qiao, J.; Wang, X.; Du, J. Impacts of changes in climate and landscape pattern on ecosystem services. Sci. Total Environ.; 2017; 579, pp. 718-728. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2016.11.036]

17. Zhang, Y.; Zhang, R.; Ma, Q.; Wang, Y.; Wang, Q.; Huang, Z.; Huang, L. A feature selection and multi-model fusion-based approach of predicting air quality. ISA Trans.; 2020; 100, pp. 210-220. [DOI: https://dx.doi.org/10.1016/j.isatra.2019.11.023]

18. Hou, L.; Wu, F.; Xie, X. The spatial characteristics and relationships between landscape pattern and ecosystem service value along an urban-rural gradient in Xi’an city, China. Ecol. Indic.; 2020; 108, 105720. [DOI: https://dx.doi.org/10.1016/j.ecolind.2019.105720]

19. Yushanjiang, Y.; Zhang, F.; Yu, H.Y.; Kung, H.T. Quantifying the spatial correlations between landscape pattern and ecosystem service value: A case study in Ebinur Lake Basin, Xinjiang, China. Ecol. Eng.; 2018; 113, pp. 94-104. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2018.02.005]

20. Shen, J.; Li, S.; Liu, L.; Liang, Z.; Wang, Y.; Wang, H.; Wu, S. Uncovering the relationships between ecosystem services and socialecological drivers at different spatial scales in the Beijing-Tianjin-Hebei region. J. Clean. Prod.; 2021; 290, 125193. [DOI: https://dx.doi.org/10.1016/j.jclepro.2020.125193]

21. Tran, D.X.; Pearson, D.; Palmer, A.; Lowry, J.; Gray, D.; Dominati, E.J. Quantifying spatial non-stationarity in the relationship between landscape structure and the provision of ecosystem services: An example in the New Zealand hill country. Sci. Total Environ.; 2022; 808, 152126. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2021.152126]

22. Foudi, S.; Spadaro, J.C.; Chiabai, A.; Polanco-Martinez, J.M.; Neumann, M.B. The climatic dependencies of urban ecosystem services from green roofs: Threshold effects and non-linearity. Ecosyst. Serv.; 2017; 24, pp. 223-233. [DOI: https://dx.doi.org/10.1016/j.ecoser.2017.03.004]

23. Peng, J.; Tian, L.; Liu, Y.X.; Zhao, M.Y.; Hu, Y.N.; Wu, J.S. Ecosystem services response to urbanization in metropolitan areas: Thresholds identification. Sci. Total Environ.; 2017; 607, pp. 706-714. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2017.06.218] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/28711000]

24. Obiang Ndong, G.; Villerd, J.; Cousin, I.; Therond, O. Using a multivariate regression tree to analyze trade-offs between ecosystem services: Application to the main cropping area in France. Sci. Total Environ.; 2021; 764, 142815. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2020.142815] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/33077232]

25. Berdugo, M.; Delgado-Baquerizo, M.; Soliveres, S.; Hernandez-Clemente, R.; Zhao, Y.C.; Gaitan, J.J.; Gross, N.; Saiz, H.; Maire, V.; Lehman, A. et al. Global ecosystem thresholds driven by aridity. Science; 2020; 367, pp. 787-790. [DOI: https://dx.doi.org/10.1126/science.aay5958] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32054762]

26. Alaboz, P.; Dengiz, O.; Demir, S.; Şenol, H. Digital mapping of soil erodibility factors based on decision tree using geostatistical approaches in terrestrial ecosystem. Catena; 2021; 207, 105634. [DOI: https://dx.doi.org/10.1016/j.catena.2021.105634]

27. Shi, G.; Zhou, Y.; Sang, Y.; Huang, H.; Zhang, J.; Meng, P.; Cai, L. Modeling the response of negative air ions to environmental factors using multiple linear regression and random forest. Ecol. Inform.; 2021; 66, 101464. [DOI: https://dx.doi.org/10.1016/j.ecoinf.2021.101464]

28. Qu, J.; Cai, W.; Zhao, Y. Learning time-dependent PDEs with a linear and nonlinear separate convolutional neural network. J. Comput. Phys.; 2022; 453, 110928. [DOI: https://dx.doi.org/10.1016/j.jcp.2021.110928]

29. An, R.; Tong, Z.; Ding, Y.; Tan, B.; Wu, Z.; Xiong, Q.; Liu, Y. Examining non-linear built environment effects on injurious traffic collisions: A gradient boosting decision tree analysis. J. Transp. Health; 2022; 24, 101296. [DOI: https://dx.doi.org/10.1016/j.jth.2021.101296]

30. Elith, J.; Leathwick, J.R.; Hastie, T. A working guide to boosted regression trees. J. Anim. Ecol.; 2008; 77, pp. 802-813. [DOI: https://dx.doi.org/10.1111/j.1365-2656.2008.01390.x]

31. Liu, M.; Chen, H.; Wei, D.; Wu, Y.; Li, C. Nonlinear relationship between urban form and street-level PM2.5 and CO based on mobile measurements and gradient boosting decision tree models. Build. Environ.; 2021; 205, 108265. [DOI: https://dx.doi.org/10.1016/j.buildenv.2021.108265]

32. Friedman, J.H. Greedy function approximation: A gradient boosting machine. Ann. Stat.; 2001; 29, pp. 1189-1232. [DOI: https://dx.doi.org/10.1214/aos/1013203451]

33. Bai, Y.; Chen, Y.; Alatalo, J.M.; Yang, Z.; Jiang, B. Scale effects on the relationships between land characteristics and ecosystem services- a case study in Taihu Lake Basin, China. Sci. Total Environ.; 2020; 716, 137083. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2020.137083] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32036149]

34. Raudsepp-Hearne, C.; Peterson, G. Scale and ecosystem services: How do observation, management, and analysis shift with scale—Lessons from Québec. Ecol. Soc.; 2016; 21, 16. [DOI: https://dx.doi.org/10.5751/ES-08605-210316]

35. Qi, Y.; Lian, X.; Wang, H.; Zhang, J.; Yang, R. Dynamic mechanism between human activities and ecosystem services: A case study of Qinghai lake watershed, China. Ecol. Indic.; 2020; 117, 106528. [DOI: https://dx.doi.org/10.1016/j.ecolind.2020.106528]

36. Mehring, M.; Ott, E.; Hummel, D. Ecosystem services supply and demand assessment: Why social-ecological dynamics matter. Ecosyst. Serv.; 2018; 30, pp. 124-125. [DOI: https://dx.doi.org/10.1016/j.ecoser.2018.02.009]

37. Chen, W.; Chi, G. Urbanization and ecosystem services: The multi-scale spatial spillover effects and spatial variations. Land Use Policy; 2022; 114, 105964. [DOI: https://dx.doi.org/10.1016/j.landusepol.2021.105964]

38. Felipe-Lucia, M.R.; Comín, F.A.; Bennett, E.M. Interactions among ecosystem services across land uses in a floodplain agroecosystem. Ecol. Soc.; 2014; 19, pp. 360-375. [DOI: https://dx.doi.org/10.5751/ES-06249-190120]

39. Roces-Díaz, J.V.; Vayreda, J.; Banqué-Casanovas, M.; Díaz-Varela, E.; Bonet, J.A.; Brotons, L.; de-Miguel, S.; Herrando, S.; Martínez-Vilalta., J. The spatial level of analysis affects the patterns of forest ecosystem services supply and their relationships. Sci. Total Environ.; 2018; 626, pp. 1270-1283. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2018.01.150]

40. Xing, L.; Zhu, Y.; Wang, J. Spatial spillover effects of urbanization on ecosystem services value in Chinese cities. Ecol. Indic.; 2021; 121, 107028. [DOI: https://dx.doi.org/10.1016/j.ecolind.2020.107028]

41. Huainan Bureau of Statistics. Huainan Statistical Yearbook 2021; China Statistics Press: Beijing, China, 2021.

42. He, T.T.; Xiao, W.; Zhao, Y.L.; Deng, X.Y.; Hu, Z.Q. Identification of waterlogging in eastern China induced by mining subsidence: A case study of Google earth engine time-series analysis applied to the Huainan coal field. Remote Sens. Environ.; 2020; 242, 111742. [DOI: https://dx.doi.org/10.1016/j.rse.2020.111742]

43. Lei, J.; Wang, S.; Wu, J.; Wang, J.; Xiong, X. Land-use configuration has significant impacts on water-related ecosystem services. Ecol. Eng.; 2021; 160, 106133. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2020.106133]

44. Tao, T.; Wang, J.; Cao, X. Exploring the non-linear associations between spatial attributes and walking distance to transit. J. Transp. Geogr.; 2020; 82, 102560. [DOI: https://dx.doi.org/10.1016/j.jtrangeo.2019.102560]

45. Huang, F.; Ochoa, C.G.; Todd Jarvis, W.; Zhong, R.; Guo, L. Evolution of landscape pattern and the association with ecosystem services in the Ili-Balkhash Basin. Environ. Monit. Assess.; 2022; 194, 171. [DOI: https://dx.doi.org/10.1007/s10661-022-09836-6] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/35147764]

46. Taylor, P.D.; Fahrig, L.; Henein, K.; Merriam, G. Connectivity is a vital element of landscape structure. Oikos; 1993; 68, pp. 571-573. [DOI: https://dx.doi.org/10.2307/3544927]

47. Amaral, Y.T.; Santos, E.M.D.; Ribeiro, M.C.; Barreto, L. Landscape structural analysis of the Lençois Maranhenses national park: Implications for conservation. J. Nat. Conserv.; 2019; 51, 125725. [DOI: https://dx.doi.org/10.1016/j.jnc.2019.125725]

48. Lee, S.W.; Hwang, S.J.; Lee, S.B.; Hwang, H.S.; Sung, H.C. Landscape ecological approach to the relationshipsof land use patterns in watersheds to water quality characteristics. Landsc. Urban Plan.; 2009; 92, pp. 80-89. [DOI: https://dx.doi.org/10.1016/j.landurbplan.2009.02.008]

49. Liu, J.; Xu, J.; Zhang, X.; Liang, Z.; Rao, K. Nonlinearity and threshold effects of landscape pattern on water quality in a rapidly urbanized headwater watershed in China. Ecol. Indic.; 2021; 124, 107389. [DOI: https://dx.doi.org/10.1016/j.ecolind.2021.107389]

50. Wang, B.; Liu, Z.; Mei, Y.; Li, W. Assessment of ecosystem service quality and its correlation with landscape patterns in Haidian District, Beijing. Int. J. Environ. Res. Public Health; 2019; 16, 1248. [DOI: https://dx.doi.org/10.3390/ijerph16071248]

51. Tu, M.; Li, W.; Orfila, O.; Li, Y.; Gruyer, D. Exploring nonlinear effects of the built environment on ridesplitting: Evidence from Chengdu. Transp. Res. Part D; 2021; 93, 102776. [DOI: https://dx.doi.org/10.1016/j.trd.2021.102776]

52. Leitão, I.A.; Ferreira, C.S.S.; Ferreira, A.J.D. Assessing long-term changes in potential ecosystem services of a peri-urbanizing Mediterranean catchment. Sci. Total Environ.; 2019; 660, pp. 993-1003. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2019.01.088]

53. Qiao, X.; Gu, Y.; Zou, C.; Xu, D.; Wang, L.; Ye, X.; Yang, Y.; Huang, X. Temporal variation and spatial scale dependency of the trade-offs and synergies among multiple ecosystem services in the Taihu Lake Basin of China. Sci. Total Environ.; 2019; 651, pp. 218-229. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2018.09.135]

54. Scholes, R.J.; Reyers, B.; Biggs, R.; Spierenburg, M.J.; Duriappah, A. Multi-scale and cross-scale assessments of social–ecological systems and their ecosystem services. Curr. Opin. Environ. Sustain.; 2013; 5, pp. 16-25. [DOI: https://dx.doi.org/10.1016/j.cosust.2013.01.004]

55. de Groot, R.S.; Alkemade, R.; Braat, L.; Hein, L.; Willemen, L. Challenges in integrating the concept of ecosystem services and values in landscape planning, management and decision making. Ecol. Complex.; 2010; 7, pp. 260-272. [DOI: https://dx.doi.org/10.1016/j.ecocom.2009.10.006]

56. McKenzie, E.; Posner, S.; Tillmann, P.; Bernhardt, J.R.; Howard, K.; Rosenthal, A. Understanding the use of ecosystem service knowledge in decision making: Lessons from international experiences of spatial planning. Environ. Plan. C Gov. Policy; 2014; 32, pp. 320-340. [DOI: https://dx.doi.org/10.1068/c12292j]

57. Chuai, X.; Huang, X.; Lai, L.; Wang, W.; Peng, J.; Zhao, R. Land use structure optimization based on carbon storage in several regional terrestrial ecosystems across China. Environ. Sci. Policy; 2013; 25, pp. 50-61. [DOI: https://dx.doi.org/10.1016/j.envsci.2012.05.005]

58. Donohue, R.J.; Roderick, M.L.; McVicar, T.R. Roots, storms and soil pores: Incorporating key ecohydrological processes into Budyko’s hydrological model. J. Hydrol.; 2012; 436–437, pp. 35-50. [DOI: https://dx.doi.org/10.1016/j.jhydrol.2012.02.033]

59. Fu, B.P. On the calculation of the evaporation from land surface (in Chinese). Chin. J. Atmos. Sci.; 1981; 5, pp. 23-31.

60. Sharp, R.; Douglass, J.; Wolny, S.; Arkema, K.; Bernhardt, J.; Bierbower, W.; Chaumont, N.; Denu, D.; Fisher, D.; Glowinski, K. et al. InVEST 3.10.2 User’s Guide 2020, The Natural Capital Project, Stanford University, University of Minnesota, The Nature Conservancy, and World Wildlife Fund. Available online: https://naturalcapitalproject.stanford.edu/software/invest (accessed on 12 November 2022).

61. Zhang, L.; Hickel, K.; Dawes, W.R.; Chiew, F.H.S.; Western, A.W.; Briggs, P.R. A rational function approach for estimating mean annual evapotranspiration. Water Resour. Res.; 2004; 40, W02502. [DOI: https://dx.doi.org/10.1029/2003WR002710]

62. Zhang, Y.; Liu, Y.; Zhang, Y.; Liu, Y.; Zhang, G. On the spatial relationship between ecosystem services and urbanization: A case study in Wuhan, China. Sci. Total Environ.; 2018; 637–638, pp. 780-790. [DOI: https://dx.doi.org/10.1016/j.scitotenv.2018.04.396]

63. Liu, W.; Zhan, J.; Zhao, C.; Zhang, F.; Teng, Y.; Chu, X.; Kumi, M.A. Spatio-temporal variations of ecosystem services and their drivers in the Pearl River Delta, China. J. Clean. Prod.; 2022; 337, 130466. [DOI: https://dx.doi.org/10.1016/j.jclepro.2022.130466]