## Materials and Methods
### Estimation of Ecosystem Services  
WR refers to the magnitude of runoff, flood, and aquifer recharge influenced by changes in ecosystems [56; 57]. We estimated WR using the Rainfall Storage Method (RSM), which considers the hydrological regulation effect of ecosystems, including forest, grassland, and cropland, compared with bare land [58; 59]. RSM can be expressed as follows:  
\[\mathrm{WR}=\mathrm{A}\times(\mathrm{J}_{0}\times\mathrm{K})\times(\mathrm{R} _{0}-\mathrm{R}_{\mathrm{g}}) \tag{7}\]  
where A is the area of the ecosystem; J\({}_{0}\) is the annual average rainfall; K is the ratio of rainfall yielding runoff to total rainfall, which is determined by precipitation and raininess; and R\({}_{0}\) and R\({}_{\mathrm{g}}\) are the runoff yield rates of bare land and the ecosystem. respectively. The thresholds of rainfall to yield runoff at different locations were collected from published literature, and the single rainfall values at the same locations were derived from three-hour TRMM precipitation data (3B42), which were revised by the measured daily precipitation data of neighboring national meteorological stations. By accumulating single rainfall values greater than the corresponding threshold, we obtained the K values at these single points. With R\({}^{2}\) = 0.79, a good linear relationship existed between these K values and annual average runoff coefficient at the same locations, thus a national spatial distribution of K value was derived from the national annual average runoff coefficient [60]. The R\({}_{\mathrm{g}}\) of forest was collected from relevant articles using meta-analysis according to factors such as vegetation types, geographic location, climate zone, and measured surface runoff. The R\({}_{\mathrm{g}}\) of grassland was calculated according to grassland coverage [61], except for alpine meadows, which adopted the study result of [PERSON] et al. [62]. The water conservation amount of floodplain wetlands and inland marshes was calculated at 8100 m\({}^{3}\) per hectare according to [PERSON] et al. [63]. SC is the reduction of soil and water loss (SWAL)--the soil erosion resulting from raindrop splash and runoff [64]--by ecosystems through their structures and processes [65; 66]. It is the difference between potential and actual soil erosion [19; 57]. Therefore, a decline in SWAL implies the improvement of SC. Using the Revised Universal Soil Loss Equation (RUSLE) [64], SWAL can be expressed as:  
\[\mathrm{SWAL}=\mathrm{R}\times\mathrm{K}\times\mathrm{LS}\times\mathrm{C} \times\mathrm{P} \tag{8}\]  
where R is the rainfall-runoff erosivity factor, K is the soil-erodibility factor, LS is the topographic factor, C is the cover-management factor, and P is the support practice factor. In this study, R was calculated using the Daily Rainfall Erosivity Model [67] based on the daily precipitation from national meteorological stations, K was obtained using the soil-erodibility nomograph [68] with the soil dataset from the Second National Soil Survey of China, C was derived according to its relationship with NDVI [69], P was obtained by reference to the USDA Handbook Number 703 and [PERSON] et al. [69], and the LS was computed using the V4.1 SRTM (Shuttle Radar Topography Mission) DEM product at a 90-m spatial resolution by the integrated method of [PERSON] et al. [70]. SP is one of the most important services supplied by ecosystems in arid and semi-arid regions, where vegetation measures help to control desertification and prevent dune movement by fixing dunes and slowing wind velocities [71; 72]. The Revised Wind Erosion Equation (RWEQ) is a combination of empirical and process-based models and has been extensively tested in simulating soil loss eroded by wind (SWIL) under broad field conditions [57; 73]. As SP is the difference between potential and actual SWIL, the less SWIL, the more SP.