Sediment Yield in Dam-Controlled Watersheds in the Pisha Sandstone Region on the Northern Loess Plateau, China

Soil erosion has become the dominant environmental issue endangering sustainable development in agriculture and the ecosystem on the Loess Plateau. Determination of watershed soil erosion rates and sediment yields is essential for reasonable utilization of water resources and soil loss control. In this study, we employed unmanned aerial vehicles (UAVs) and structure-from-motion (SfM) photogrammetry to determine the sediment yields in 24 dam-controlled watersheds in the Pisha sandstone region of the northern Loess Plateau. High differences in total sediment were trapped before the check dams due to their running periods and sediment yields. The estimated specific sediment yield ranged from 34.32 t/(ha·a) to 123.80 t/(ha·a) with an average of 63.55 t/(ha·a), which indicated that the Pisha sandstone region had an intense soil erosion rate. Furthermore, the modified Sediment Distributed Delivery (SEDD) model was applied to identify the erosion-prone areas in the watersheds, and the sediment retained in the check dams were used for model calibration. The performance of the model was acceptable, and the modeling results indicated that the steep Pisha sandstone was the major sediment source for the watersheds, accounting for approximately 87.37% of the sediment yield. Catchment area, erosive precipitation, and badland proportion were the key factors for sediment yield in the dam-controlled watersheds of the Pisha sandstone region, according to multiple regression analyses. These findings indicated that the modified SEDD model is very efficient in identifying spatial heterogeneities of sediment yield in the watershed but requires comprehensive calibration and validation with long-term observations. The Pisha sandstone region is still the key area of soil erosion control in the Loess Plateau, which needs more attention for soil and water conservation due to high sediment yield.


Introduction
Soil erosion is a significant environmental issue that undermines the sustainable development in agriculture and ecosystems of the Loess Plateau [1,2]. Soil erosion results in loss of fertile surface topsoil, leading to degradation of soil quality and reduction of the limited amount of available arable land, thereby undermining food and environmental security [3,4]. Additionally, a considerable amount of sediment from soil erosion could cause serious sedimentation in reservoirs or riverbeds, thus affecting the regulation of water resources and the security of water transportation [2,4].
Check dams, as one of the most important gully control engineering works, can stabilize slopes, mitigate soil erosion, and reduce sediment transportation downstream [5][6][7]. Numerous studies have paid close attention to the impacts of check dams on the hydrological and geomorphological processes in the past decades [8][9][10]. Remaître et al. [11] proposed that the debris-flow intensity might decrease efficiency by the number of check dams near the source area. However, some studies found that check-dams might exacerbate local erosion and affect the watershed sediment budget [12,13]. Boix-Fayos et al. [14] proposed that the impact of check dams was negative on soil erosion mitigation because the sedimentary dynamics and riverbed stability were altered by the check dams.
Additionally, many studies address that sedimentation before the check dams could provide an effective approach to estimating sediment yield in watersheds [15,16]. Li et al. [17] estimated the relationship between sediment yield and storm characteristics by digging sedimentation profiles behind four dams on the Loess Plateau. Porto et al. [18] suggested that 137 Cs and 210 Pb ex measurements could offer an effective method for estimating both erosion and sediment redistribution in a small forested watershed. Recently, field surveys were employed in combination with high-precision differential GPS techniques to estimate the sediment yield of the small watersheds [19]. Alfonso-Torreno et al. [20] applied unmanned aerial vehicles (UAVs) and structure-from-motion (SfM) photogrammetry together with a geometric method to estimate sediment volumes trapped by check dams and demonstrated satisfactory agreement between estimated and measured sediment volumes. These previous studies confirmed that sedimentation in the check dams provided important information to estimate the overall sediment yields in watersheds.
Spatial patterns of soil erosion and sediment yield are commonly obtained from soil erosion modeling [21][22][23]. Liu et al. [24] applied the SWAT model to determine the watershed sediment source and found that channel erosion contributed about 60% of the total sediment in the watershed. Boakye et al. [25] indicated 21.3% of the watershed was an erosion-prone area in Pra River Basin using the SEDD model combined with GIS. These studies indicated that soil erosion models were useful tools to estimate spatial heterogeneities in soil erosion and could be applied to determine soil erosion-prone areas in a macroscale catchment.
Considerable studies have addressed a significant sediment load reduction in the middle and lower reaches of the Yellow River due to the implementation of numerous soil and water conservation measures over the past 60 years [26,27]. Yao et al. [28] also found sediment reduction up to 60% in the Hekouzhen-Longmen section from 2000 to 2012 due to the implementation of soil and water conservation measures. The Pisha sandstone region is the most seriously erosion-prone area on the Loess Plateau, contributing 71.1% of the coarse sediment into the reaches of the Yellow River [29,30]. By now, spatial patterns of soil erosion in this region have not been comprehensively investigated. Sediment yields are not clear due to very few studies in the Pisha sandstone area.
Thus, the aims of this study were to (i) quantify the sediment yields in dam-controlled watersheds in the Pisha sandstone region based on UAVs and SfM; (ii) identify the key soil erosion-prone area using the calibrated SEDD model; and (iii) obtain the dominant environment factors related to sediment yield.

Study Area
The study was conducted in the Pisha sandstone area in the northern Loess Plateau of China (38 • 10 -40 • 10 N, 108 • 45 -111 • 31 E), which is mainly located in the border area of Shanxi, Shaanxi, and Mongolia. The Pisha sandstone area covers an area of 1.67 × 10 4 km 2 and is dominated by wind and water erosion ( Figure 1). The climate is temperate continental arid and semi-arid monsoon climate with a mean annual precipitation of 280-400 mm. The annual rainfall is mainly concentrated from June to September with high-intensity storms, which decrease from southeast to northwest in terms of spatial variation. The landscape is mainly in the form of dense gullies with sparse vegetation. The main soil types in this region are loess, desert sand, and coarse weathered sandstone-locally called Pisha sandstone (Figure 2a), which is characterized by loose structure, poorly bonded mechanisms, weak antiscourability, and erosion resistance [31]. Numerous soil and water conservation measures

Selection of Watersheds
In our case, we selected the check dam-controlled watersheds using the Google Earth images together with a field survey. The selection of check dams required the following conditions: (i) Check dams were located upstream of the channels, ensuring accurate estimation of sediment trapped behind check dams. (ii) Check dams had a certain silt period, and sedimentation cores could be easily obtained to estimate sediment yield. (ⅲ) The selected check dams had no spillway or drainage structure. In total, 24 dam-controlled watersheds were selected to measure sediment retained by the check dams (Table 1).  Pisha sandstone (Figure 2a), which is characterized by loose structure, poorly bonded mechanisms, weak antiscourability, and erosion resistance [31]. Numerous soil and water conservation measures have been carried out in this region, where check dams have become the dominant engineering measures (Figure 2b).

Selection of Watersheds
In our case, we selected the check dam-controlled watersheds using the Google Earth images together with a field survey. The selection of check dams required the following conditions: (i) Check dams were located upstream of the channels, ensuring accurate estimation of sediment trapped behind check dams. (ii) Check dams had a certain silt period, and sedimentation cores could be easily obtained to estimate sediment yield. (ⅲ) The selected check dams had no spillway or drainage structure. In total, 24 dam-controlled watersheds were selected to measure sediment retained by the check dams (Table 1).

Selection of Watersheds
In our case, we selected the check dam-controlled watersheds using the Google Earth images together with a field survey. The selection of check dams required the following conditions: (i) Check dams were located upstream of the channels, ensuring accurate estimation of sediment trapped behind check dams. (ii) Check dams had a certain silt period, and sedimentation cores could be easily obtained to estimate sediment yield. (iii) The selected check dams had no spillway or drainage structure. In total, 24 damcontrolled watersheds were selected to measure sediment retained by the check dams (Table 1).

Sediment Yield Estimations for Small Watersheds
Sediment yield can be estimated by sediment volume and bulk density of sedimentation. As shown in Figure 3, we took images of the watersheds by an unmanned air vehicle (UAV, DJI Phantom 4 professional, https://www.dji.com, (accessed on 1 October 2019)) at an approximate altitude of 150 m. UAV images were processed by Agisoft PhotoScan Professional software (v. 1.4.5) to generate orthophotographs ( Figure 3b) and DEM ( Figure  3c) and were georeferenced by ground control points with RTK. DEMs were imported into ArcGIS 10.5 (http://www.esri.com, (accessed on 15 October 2019)) for hydrological analysis to obtain the watershed boundary, and the extent of the sediment wedge was delineated according to the orthophotographs. At the same time, we used manual drilling to obtain the sedimentation depth ( Figure 3a). Due to the steep slope of the selected dam-controlled watersheds, the sediment wedge was regarded as an irregular prism (the deep siltation area is almost the same as the extent of the wedge surface). Sediment volume for each watershed was calculated by the raster calculator function in ArcGIS 10.5. Three samples were collected from the sedimentation profiles in each dam-controlled watershed by using steel rings (100 cm 3 ) to estimate the mean dry bulk density (Table 1). Since the selected check dams do not have a spillway or sluicing gate, the incoming upstream sediment has been completely deposited before the check dams. Thus, the annual sediment yield (SY, t/a) for the dam-controlled watersheds was estimated by using the following Equation (1): where V is the sediment volume of the sediment wedge (m 3 ), ρ is the mean bulk density of the sedimentation (kg/m 3 ), T is the silting period (a). The specific sediment yield (SSY, t/(ha•a)) was then estimated by: where A is the drainage area (ha).

The SEDD Model
The SEDD model, as a spatially distributed approach, was applied to estimate sediment yield in the watershed. It coupled the RUSLE model [32] and the sediment delivery ratio (SDR) for each watershed morphological unit [33,34]. The mean annual gross soil erosion (t/(ha•a)) for each grid cell was calculated by: where Ri is the rainfall erosivity factor (MJ mm/ha/h); Ki is the soil erodibility factor (t ha h/(MJ ha mm)); Li and Si are topographic factors, which are dimensionless variables; Ci represents the cover management factor, and Pi represents support practice factors, which are dimensionless variables. Due to the absence of 30-min rainfall intensity data, the method proposed by Zhang et al. [35] was applied to obtain the R values based on the daily erosive precipitation (>12 mm) from 11 rainfall stations around the dam-controlled watershed from 2006 to 2018 ( Figure 1). Since the selected check dams do not have a spillway or sluicing gate, the incoming upstream sediment has been completely deposited before the check dams. Thus, the annual sediment yield (SY, t/a) for the dam-controlled watersheds was estimated by using the following Equation (1): where V is the sediment volume of the sediment wedge (m 3 ), ρ is the mean bulk density of the sedimentation (kg/m 3 ), T is the silting period (a). The specific sediment yield (SSY, t/(ha·a)) was then estimated by: where A is the drainage area (ha).

The SEDD Model
The SEDD model, as a spatially distributed approach, was applied to estimate sediment yield in the watershed. It coupled the RUSLE model [32] and the sediment delivery ratio (SDR) for each watershed morphological unit [33,34]. The mean annual gross soil erosion (t/(ha·a)) for each grid cell was calculated by: where R i is the rainfall erosivity factor (MJ mm/ha/h); K i is the soil erodibility factor (t ha h/(MJ ha mm)); L i and S i are topographic factors, which are dimensionless variables; C i represents the cover management factor, and P i represents support practice factors, which are dimensionless variables. Due to the absence of 30-min rainfall intensity data, the method proposed by Zhang et al. [35] was applied to obtain the R values based on the daily erosive precipitation (>12 mm) from 11 rainfall stations around the dam-controlled watershed from 2006 to 2018 ( Figure 1). The K factor was obtained through field sampling at 10-15 sites for each watershed. The K values were calculated based on the erodibility nomograph method [36], which was estimated according to soil properties such as soil texture, organic matter content, and permeability.
The study area was characterized by complex topography and geomorphology, and large areas were covered by steep gullies with gradients higher than 20%. Thus, the method proposed by Liu et al. [37] was applied to calculate the LS factor, modified from the RUSLE model.
The C factor demonstrates integrated influences of cropping and management practices on soil erosion [38]. Low C factor values represent areas with dense vegetation, such as forestland, and higher values associated with bare land without management practices. We assigned C factor values to different land-use types according to a combination of field-based expertise and the existing literature (Table 2) [39,40]. The P-factor reflects the impact of practices such as contouring, strip-cropping, terraces, or contour furrows on soil erosion. In this study, field observation showed that the only water and soil conservation measure for each watershed was fish-scale pits. Therefore, we assigned a P-factor value of 0.187 to forestland, and the P-factor of other land-use types was 1.
The SEDD model divides a watershed into geomorphic units, and the SDR for each unit is calculated. The SDR i was estimated by the following equation: where β is a constant coefficient with no dimension for the specific watershed and t i is the movement time (h) for each cell. M is the number of cells located in the flow paths, l i is the flow length (m), and v i is the stream speed (m/s), which is estimated by: where s i is the slope of the cell, and d i is a parameter closely linked to land use (m/s). The d i values can be obtained from Ferro and Porto [34] and Ferro et al. [41]. The constantcoefficient β depends critically on the basin morphology and can be estimated using an inverse modeling algorithm. Fernandez et al. [42] proposed that the SDR w could be calculated according to the relationships between sediment yield and drainage area: where S w is the catchment area (km 2 ), and b is a dimensionless coefficient that usually takes a value of 0.0328 for RUSLE [41]. When SDR w is determined, β can be calculated by Equation (7): Land 2021, 10, 1264 7 of 17 where N is the total amount of cells for the watershed, l i and s i are the flow length and slope for each cell, respectively, and a i is the cell area. The sediment yield (t/a) for each watershed is calculated by: where N represents the total amount of cells for the watershed, SDR i is the SDR for each cell, and A i represents the soil erosion rate calculated by Equation (4). Land use/cover, topography, climatic and soil type data were required to run the model. Table 3 lists the data sources. All input data were pre-processed into spatially consistent raster layers to run the SEDD model. We used a spatial resolution of 2 × 2 m, corresponding to the resampled DEM for each dam-controlled watershed, which was derived from the UAVs and SfM. The model was calibrated by comparing the modeling results with observations from field measurements. We used the Nash-Sutcliffe efficiency (NSE) values, and relative root mean squared error (RRMSE) to evaluate the model performance. These indices are expressed as follows: where n represents the number of observations, O i represents the measured sediment yield for each watershed, O mean represents the average measured values, and M i represents the modeled values for each watershed.

Statistical Analysis
The Pearson correlation coefficients were used to explore the relationships between the sediment yield and environmental variables (Table 4). Correlation coefficients may reveal only a partial relationship between sediment yield and environmental variables. Therefore, a multivariate regression analysis was applied to estimate sediment yield with a number of variables. Before establishing the final multivariate regression model, we calculated the variance inflation factor (VIF) to reduce the collinearity variables, which followed the "rule of 10" [43]. If the value of VIF is greater than 10, the corresponding variable was removed from the model. Finally, we established the optimal multivariate regression model to explain the sediment yield.  Figure 4 shows the estimated annual sediment yield (SY) and specific sediment yield (SSY) for each dam-controlled watershed based on the UAVs and SfM. We found that a significant difference existed in the total sediment for the selected watersheds due to their running periods and land surface characteristics. The SYs of 24 dam-controlled watersheds varied from 82.45 t/a to 8223.98 t/a, with an average of 1948.16 t/a (Figure 4a). Additionally, large discrepancies existed among the SSY in the dam-controlled watersheds. Watershed SSYs ranged from 34.32 t/(ha·a) to 123.80 t/(ha·a) (Figure 4b). The average specific sediment yield of the selected watersheds is 63.55 t/(ha·a), which indicated that the Pisha sandstone region had an intense soil erosion rate. These findings indicate that the Pisha sandstone area is still the key point of soil erosion control in the Loess Plateau due to high specific sediment yields.

Model Calibration and Validation
The SEDD model was calibrated according to the sedimentation in the check dams. The values of NSE and RRMSE were used to evaluate the model performance. Figure 5

Model Calibration and Validation
The SEDD model was calibrated according to the sedimentation in the check dams. The values of NSE and RRMSE were used to evaluate the model performance. Figure 5 shows the relationship between the measured and simulated specific sediment yield (SSY) for each dam-controlled watershed. In general, the model provides relatively satisfactory results (NSE = 0.85, RRMSE = 0.83). The simulated SSY varied from 24.79 to 131.85 t/(ha·a), and the measured SSY varied from 34.32 t/(ha·a) to 123.80 t/(ha·a). However, the model was inclined to underestimate the SSY compared with the measured values. The average values of modeled and measured SSY are 59.84 t/(ha·a) and 63.55 t/(ha·a), respectively.

Model Calibration and Validation
The SEDD model was calibrated according to the sedimentation in the check The values of NSE and RRMSE were used to evaluate the model performance. Fi shows the relationship between the measured and simulated specific sediment yield for each dam-controlled watershed. In general, the model provides relatively satisf results (NSE = 0.85, RRMSE = 0.83). The simulated SSY varied from 24.79 to 131.85 t/ and the measured SSY varied from 34.32 t/(ha•a) to 123.80 t/(ha•a). However, the was inclined to underestimate the SSY compared with the measured values. The av values of modeled and measured SSY are 59.84 t/(ha•a) and 63.55 t/(ha•a), respectiv  Figure 6 shows sediment yield from different land use/cover types for the se watersheds. On average, the steep bare Pisha sandstone accounted for 21.28% of th of the watershed (Figure 6a) but contributed approximately 87.37% of the sedime the remaining 12.63% of the sediment originated from grassland, shrubland, and plain. Besides, shrubland and fluvial plain accounted for only a small part of the sed (2.77% and 2.10%, respectively.) (Figure 6b). These findings suggest that the sedim tion in the check dams mainly originated from the steep gullies with weathered sand Figure 5. Comparison of the measured and simulated specific sediment yields. Figure 6 shows sediment yield from different land use/cover types for the selected watersheds. On average, the steep bare Pisha sandstone accounted for 21.28% of the area of the watershed (Figure 6a) but contributed approximately 87.37% of the sediment and the remaining 12.63% of the sediment originated from grassland, shrubland, and fluvial plain. Besides, shrubland and fluvial plain accounted for only a small part of the sediment (2.77% and 2.10%, respectively.) (Figure 6b). These findings suggest that the sedimentation in the check dams mainly originated from the steep gullies with weathered sandstone and that the steep Pisha sandstone region was the hotspot of erosion and sediment yield for each watershed. and that the steep Pisha sandstone region was the hotspot of erosion and sediment yield for each watershed. To better identify the spatial distribution of sediment yield, we compared the specific sediment yield of two dam-controlled watersheds (ID = 8 and ID = 21). The land-use types in the two watersheds are mainly grassland, shrubland, fluvial plain, and sandstone (Figure 7a,b), whereas we can clearly distinguish the discrepancy in sediment yield between watersheds 8 and 21 (Figure 7c,d). High sediment yields (>100 t/(ha•a)) are mainly distributed in steep gullies, and the low values (<25 t/(ha•a)) are mainly located in the vegetated regions with gentle slopes. Additionally, the spatial distribution of soil erosion demonstrates that bare, weathered sandstone in steep gullies is the most erosion-prone area in the watershed, producing the highest sediment yield in both watersheds. To better identify the spatial distribution of sediment yield, we compared the specific sediment yield of two dam-controlled watersheds (ID = 8 and ID = 21). The land-use types in the two watersheds are mainly grassland, shrubland, fluvial plain, and sandstone (Figure 7a,b), whereas we can clearly distinguish the discrepancy in sediment yield between watersheds 8 and 21 (Figure 7c,d). High sediment yields (>100 t/(ha·a)) are mainly distributed in steep gullies, and the low values (<25 t/(ha·a)) are mainly located in the vegetated regions with gentle slopes. Additionally, the spatial distribution of soil erosion demonstrates that bare, weathered sandstone in steep gullies is the most erosion-prone area in the watershed, producing the highest sediment yield in both watersheds.    Figure 8 shows the Pearson's correlation coefficient (r) among all variables. Circles of different sizes and colors represent the degree and sign of correlation. SY in dam-controlled watersheds was positively correlated with CA (0.85), SPI (0.78), TR (0.67), SLO (0.61), and BL (0.60), which are significant at a confidence level of 99%. Furthermore, SY was correlated with GC (0.43), ProC (0.42), R (0.46), and VC (0.51) (p < 0.05). However, we found collinearity among environmental variables such as GC and ProC (0.98), TR and SLO (0.96), TR and BL (0.84), SLO and BL (0.8), etc. This evidence indicated that the correlation between environmental variables and SY could be enhanced or offset by the correlation among these explanatory variables. Before establishing the final multivariate regression model, the variance inflation factor (VIF) was used to reduce some collinearity variables. When the VIF values of environment variables are greater than 10, the environment variables will be excluded from the multiple regression analysis for SY. Thus, ProC (VIF = 100.42), SLO (VIF = 21.47), and TR (VIF = 10.20) were removed from the multiple regression analysis for SY. Finally, we established a multiple regression model for sediment yield (Equation (11)). The regression results showed that the CA, R, and BL were positively correlated with SY. The equation also indicated that catchment area, erosive precipitation, and badland proportion were significant driving factors for sediment yield in the small watersheds of the Pisha sandstone region (Adjusted R 2 = 0.90). SY = 0.72CA + 0.16R + 0.35BL (n = 24, R 2 = 0.90)

Estimation of Sediment Yield in The Small Watersheds
In our case, we quantified sediment yield in ungauged dam-controlled watersheds Before establishing the final multivariate regression model, the variance inflation factor (VIF) was used to reduce some collinearity variables. When the VIF values of environment variables are greater than 10, the environment variables will be excluded from the multiple regression analysis for SY. Thus, ProC (VIF = 100.42), SLO (VIF = 21.47), and TR (VIF = 10.20) were removed from the multiple regression analysis for SY. Finally, we established a multiple regression model for sediment yield (Equation (11)). The regression results showed that the CA, R, and BL were positively correlated with SY. The equation also indicated that catchment area, erosive precipitation, and badland proportion were significant driving factors for sediment yield in the small watersheds of the Pisha sandstone region (Adjusted R 2 = 0.90). SY = 0.72CA + 0.16R + 0.35BL (n = 24, R 2 = 0.90) (11) Land 2021, 10, 1264 12 of 17

Estimation of Sediment Yield in The Small Watersheds
In our case, we quantified sediment yield in ungauged dam-controlled watersheds based on the UAVs and SfM. The proposed method could allow rapid acquisition of DEM and orthophotos at different scales, overcoming inaccurate sediment estimation due to the insufficient quality of coarse resolution images generated by other measurement techniques. This method could reduce numerous field studies compared to sediment estimation via long-term dynamic monitoring or traditional ground-based measurements. However, the major limitation of this method was the inaccuracy of the SfM techniques to reconstruct topography in areas with dense vegetation cover due to the occlusion caused by vegetation [44,45]. Alfonso-Torreno et al. [20] proposed that the traditional geomorphological techniques, such as Total Stations and GNSS, could be used as an alternative to reconstruct the topography of dam-controlled watersheds with dense vegetation cover. Besides, we believed another main uncertainty source for sediment estimation was that the sediment wedge was generalized as an irregular prism. Due to the complicated geomorphological characteristics of the gully before the dam, sedimentation volume might be over/underestimated due to the great variability in shape and thickness of the sediment wedge. Future works should incorporate advanced geological exploration techniques to improve sediment volume estimation.
As mentioned above, a large variability existed in estimated specific sediment yield for each dam-controlled watershed. However, the range of estimated specific sediment yield was close to previous studies on the Loess Plateau (Table 5). We compared our estimation with some related studies to better understand the differences in specific sediment yield for small watersheds. Overall, the average specific sediment yield was slightly lower than that of previous studies in the Pisha sandstone area [46][47][48][49]. This discrepancy might be attributable to different land use/cover and large-scale afforestation projects. Fu et al. [2] suggested soil erosion and sediment yield reduced greatly due to the impact of the "Grainfor-Green" project. However, the average specific sediment yield was higher than that in the loess region on the Loess Plateau [31,50]. The differences between our estimation and their results may be caused by diverse landscapes and the widely distributed bare, weathered sandstone in the study areas. This was also confirmed by Jiao et al. [51] in the northern Loess Plateau.

Model Performance
The traditional empirical model (RUSLE) does not specifically account for gully erosion or steep slope zones. The Loess Plateau is characterized by complicated topography with dense gullies, where high-intensity storms may result in the frequent occurrence of landslides and other disasters. Gully erosion is one of the important soil erosion processes on the Loess Plateau, causing a large amount of sediment into the river. However, the RUSLE model can not take into account these soil erosion processes due to its limitations [52]. The SEDD model could not only identify the spatial distribution of sediment yield and determine the watershed sediment source but also required less data and time to run due to its easy parameterization. Therefore, we applied the SEDD model to identify the erosion-prone area for the selected watersheds.
In our case, the model was calibrated according to the sedimentation in the check dams. The model performed satisfactorily with high model efficiencies, but most simulated values were slightly lower than those measured from the field survey. This difference might be associated with inaccurate estimations of sediment yield for small watersheds due to complex erosion processes or errors from sediment volume estimation. Bussi et al. [53] found the watershed sediment yield estimates based on reservoir sedimentation were often questionable for the reservoir capacity estimation without initial topography. De Vente et al. [54] argued that errors such as volume calculations might lead to an overestimation or underestimation of the watershed sediment yield. Much work should be done in the future to improve model performance; for example, a large number of field measurements are implemented to verify the input parameters of the model because the quality of the simulation results greatly depend on these parameters. Besides, by supplementing the sampling frequency of check dams and construction of sediment time series, the model performance will be enhanced by a comprehensive validation with long-term observations.

The Key Factors Affecting Sediment Yield
The shape and geomorphological features of the watershed affected sediment yield mainly by influencing energy transformation and sediment movement [55]. This study evidenced a significant positive relationship between the catchment area and sediment yield. The explanations for this positive correlation are as follows: (i) high topographic relief in small watersheds led to shorter distances between sediment sources and channels, resulting in high SDR in the watershed; (ii) soil erosion and sediment yield processes in the watershed were accelerated by localized storms. Furthermore, the sediment transportation would be improved due to the distribution of considerable steep gullies linking the hillslope and the riverbed in the watershed.
Rainfall amount and intensity are important characteristics of precipitation that affect runoff generation and soil detachment, which made an impact on erosion process in the watershed [56,57]. Additionally, soil erosion was triggered by erosive precipitation, which referred to precipitation events that exceed the threshold for sediment-producing precipitation [58]. In our case, rainfall erosivity (R) was used as a comprehensive indicator to characterize the effects of erosive precipitation on sediment yield. The regression result indicated a significant positive relationship between R and sediment yield, which can be easy to understand because sediment produced by soil erosion requires runoff [59,60].
The badland, as a major sediment source, was found to also be the dominant factor that influenced the watershed sediment yield [61][62][63][64]. We also found a positive correlation between badland proportion and sediment yield. This positive correlation was attributed to the widespread distribution of bare, weathered sandstone. Weathered sandstone is characterized by weak cementation ability [65], poor anti-scourability [66], and high erodibility. Therefore, it can be disintegrated rapidly by water, leading to frequent gully erosion and considerable sediment yield [67,68]. Besides, high-intensity human activities, such as overgrazing or indiscriminate mining, majorly damaged the ecological environment, and the badland proportion increased sharply due to poor vegetation cover, thus resulting in severe erosion and high sediment yields in the Pisha sandstone region [69].

Conclusions
In our case, we applied the UAVs and SfM and the modified SEDD model for estimating the sediment yield in 24 dam-controlled watersheds in the Pisha sandstone region on the northern Loess Plateau and discussed the dominant factors on sediment yield in the watershed. The main conclusions are as follows.
Large differences in total sediment were retained in the check dams due to different running periods and sediment yields. The specific sediment yields of the watersheds varied from 34.32 t/(ha·a) to 123.80 t/(ha·a) with an average of 63.55 t/(ha·a), indicating that the Pisha sandstone region had an intense soil erosion rate. The SEDD model was verified through the sedimentation in the check dams. A good agreement was found between the modeling values and observations, which demonstrated that the model performance was acceptable and could be used to identify the erosion-prone area in the watershed. Additionally, approximately 87.37% of the sediment was contributed by the weathered sandstone in the steep gullies.
The final multiple regression model for sediment yield (SY) included catchment area (CA), rainfall erosivity factor (R), and badland proportion (BL) were positively associated with SY. The results indicated that catchment area, erosive precipitation, and badland proportion were the dominant factors for sediment yield in the dam-controlled watersheds of the Pisha sandstone region.
Based on our study, the modified SEDD model is very efficient for estimating the spatial distribution of sediment yield but requires comprehensive calibration and validation. The Pisha sandstone region needs comprehensive soil and water conservation measures due to its high sediment yields in the Loess Plateau.
Author Contributions: X.M. and G.Z. conceptualized this study; F.X. and G.Z. designed the study; P.G. and W.S. performed data checks; F.X. and G.Z. wrote the original draft preparation; G.Z. and P.T. edited and revised the original manuscript. All authors have read and agreed to the published version of the manuscript.