Sensitivity Assessment of Spatial Resolution Difference in DEM for Soil Erosion Estimation Based on UAV Observations: An Experiment on Agriculture Terraces in the Middle Hill of Nepal

: Soil erosion in the agricultural area of a hill slope is a fundamental issue for crop productivity and environmental sustainability. Building terrace is a very popular way to control soil erosion, and accurate assessment of the soil erosion rate is important for sustainable agriculture and environmental management. Currently, many soil erosion estimations are mainly based on the freely available medium or coarse resolution digital elevation model (DEM) data that neglect micro topographic modiﬁcation of the agriculture terraces. The development of unmanned aerial vehicle (UAV) technology enables the development of high-resolution (centimeter level) DEM to present accurate topographic features. To demonstrate the sensitivity of soil erosion estimates to DEM resolution at this high-resolution level, this study tries to evaluate soil erosion estimation in the Middle Hill agriculture terraces in Nepal based on UAV derived high-resolution (5 × 5 cm) DEM data and make a comparative study for the estimates by using the DEM data aggregated into different spatial resolutions (5 × 5 cm to 10 × 10 m). Firstly, slope gradient, slope length, and topographic factors were calculated at different resolutions. Then, the revised universal soil loss estimation (RUSLE) model was applied to estimate soil erosion rates with the derived LS factor at different resolutions. The results indicated that there was higher change rate in slope gradient, slope length, LS factor, and soil erosion rate when using DEM data with resolution from 5 × 5 cm to 2 × 2 m than using coarser DEM data. A power trend line was effectively used to present the relationship between soil erosion rate and DEM resolution. The ﬁndings indicated that soil erosion estimates are highly sensitive to DEM resolution (from 5 × 5 cm to 2 × 2 m), and the changes become relatively stable from 2 × 2 m. The use of DEM data with pixel size larger than 2 × 2 m cannot detect the micro topography. With the insights about the inﬂuencing mechanism of DEM resolution on soil erosion estimates, this study provides important suggestions for appropriate DEM data selection that should be investigated ﬁrst for accurate soil erosion estimation.


Introduction
Land is one of the most important natural resources for agricultural activities to support people's livelihood [1]. As being a crucial issue of decreasing productivity, the loss of fertile soil induced by soil erosion frequently happens in hilly slopes [2], and the effects of soil erosion go beyond the loss of fertile land [3]. The global soil erosion rate was estimated at 35.9 Pg yr −1 in 2012, and there is a potentially increasing trend in global soil erosion driven by cropland expansion [4]. Under the background of climate change, soil ISPRS Int. J. Geo-Inf. 2021, 10, 28 2 of 17 erosion has been and will always be a key issue for the sustainable development at the framework of future earth [3,5]. Meanwhile, the catastrophic floods that have frequently happened in recent years greatly induced big amount of soil erosion, posing a big pressure on local agricultural development [6,7].
Usually, soil erosion process by water in the hill slope is highly controlled by its terrain features [8]. Therefore, increasing soil water retention capacity by improving infiltration and reducing peak discharge rate of surface runoff has been well recognized as the main direction for soil and water conservation. Terracing cultivated land in the hill slope can obviously achieved this effect [9]. Therefore, terrace is one of the most predominant forms of agriculture in the Asia and the Pacific regions [10].
As a typical mountain country in the South Asia, Nepal has experienced serious degradation of soil resources over the past three decades especially in the Middle Hills as a result of the expansion of agricultural land and the increase in cropping intensity [2]. The national land degradation assessment report has revealed that 3.16 million hectares (11.81% of the total area of Nepal) have been affected by the degradation process. To preserve soil and water resources, both irrigated contour terraces and rain-fed outward sloping terraces are found in most of the hill slopes of Nepal. The modification of topographic features at small scales by these terraces greatly controls soil erosion [11], which is a long tradition of soil erosion controlling technique in the cultivated hill slope of Nepal. Farmers maintain their agriculture terraces regularly. They also grow trees at rain-fed agriculture terraces to stabilize it and fodder supply for their livestock [12]. These are the major characteristics of integrated subsistence farming system of agriculture hill slopes of Nepal.
From field measurements, it has been confirmed that the controlling effect of soil erosion is quite different for different types of terraces [13]. Thus, accurate detection of real topographic features of hilly slopes is crucial for soil erosion monitoring and regulation in the slope. With the development of digital elevation model (DEM), data are a powerful vehicle for conveying essential surface topographic information, which is very useful for terrain analysis of physical surface including various hydrological and biophysical properties [14]. Therefore, the resolution of the DEM data plays an important role in controlling the soil erosion assessment. Currently, many models have been designed for soil erosion estimation, such as the universal soil loss equation (USLE), modified universal soil loss equation (MUSLE), and revised universal soil loss equation (RUSLE). The RUSLE model estimates sheet and rill erosion caused by raindrop impact [15], and this model does not cover the gully erosion. With the combination of geographic information system (GIS) and remote sensing technologies, RUSLE estimates soil erosion loss on a cell-by-cell basis. Soil and water assessment that uses DEMs with the resolutions from 5 × 5 to 14 × 14 m reveals the decreasing trend when DEMs were resampled to coarser [16]. Similar results were achieved by different studies with higher accuracies of higher resolution DEM as compared to coarse resolution [17,18]. Thus, there are the significant effects of spatial resolution of DEM on topography-based geomorphological and hydrological assessment [19][20][21][22][23].
From current studies, it can be concluded that most of them were conducted using the existing open-source DEMs, which range from 10 × 10 m to kilometer level. There is a lack of related discussions or analyses for high spatial resolution (sub-meter level) due to the limitations of such high-resolution DEM data. In recent years, the development of unmanned aerial vehicle (UAV) technique provides an easy tool to derive high resolution DEM with the spatial resolution of centimeter level [24]. Thus, the UAV-based high spatial resolution DEM offers valuable information of surface topography especially land surface processes and hydrological analysis. UAV-acquired high-resolution DEM has been popular for several geomorphological and hydrological analyses in recent decades [25][26][27], because it provides an important opportunity to estimate soil erosion at this high-resolution level. Therefore, with a view to accurate estimation of the soil erosion in Nepal, this study has been conducted to investigate the sensitivity of the estimated of soil erosion rates based on the RUSLE model for agriculture terraces in the hill slopes of Nepal. The findings of the study will be helpful to understand how far the spatial resolution of DEM is sensitive for the estimation of soil erosion rates and the values of topographic factors for geomorphic and hydrologic processes research, especially in the narrow agriculture terraces in the hill slopes such as those in Nepal.

Study Area
The study area is located in the Marshynangdi river basin that lies in the hilly region of the central part of Nepal. The extension of the study area is 28 • 12 33 north to 28 • 13 07 north of latitude and 84 • 25 28 east to 84 • 26 37 east of longitude ( Figure 1) with the elevation range from 1072 to 1593 m from the mean sea level. It has a warm temperate climate with an annual average rainfall of 3374 mm (more than 80% fall in summer season from May to August). Therefore, the summer season usually experiences the most serious soil erosion in this region. resolution level. Therefore, with a view to accurate estimation of the soil erosion in Nepal, this study has been conducted to investigate the sensitivity of the estimated of soil erosion rates based on the RUSLE model for agriculture terraces in the hill slopes of Nepal. The findings of the study will be helpful to understand how far the spatial resolution of DEM is sensitive for the estimation of soil erosion rates and the values of topographic factors for geomorphic and hydrologic processes research, especially in the narrow agriculture terraces in the hill slopes such as those in Nepal.

Study Area
The study area is located in the Marshynangdi river basin that lies in the hilly region of the central part of Nepal. The extension of the study area is 28°12'33″ north to 28°13'07″ north of latitude and 84°25'28″ east to 84°26'37″ east of longitude ( Figure 1) with the elevation range from 1072 to 1593 m from the mean sea level. It has a warm temperate climate with an annual average rainfall of 3374 mm (more than 80% fall in summer season from May to August). Therefore, the summer season usually experiences the most serious soil erosion in this region. Farmers manage their cultivated land by making complete horizontal terraces at irrigated plots and outward sloping terraces for rain-fed plots. They grow grasses and fodder trees at terrace riser. Intensive subsistence farming practice is evident in the area as intensive use of cultivated land, which integrates livestock, forest, and grasses in their system of livelihood. As shown in Figure 1, six plots were selected as the study areas. Among them, three were rain-fed and three were irrigated agriculture plots. In local term, Farmers manage their cultivated land by making complete horizontal terraces at irrigated plots and outward sloping terraces for rain-fed plots. They grow grasses and fodder trees at terrace riser. Intensive subsistence farming practice is evident in the area as intensive use of cultivated land, which integrates livestock, forest, and grasses in their system of livelihood. As shown in Figure 1, six plots were selected as the study areas. Among them, three were rain-fed and three were irrigated agriculture plots. In local term, a rain-fed cultivated terrace is known as bari, and an irrigated terrace is known as khet. Therefore, the selected rain-fed plots were named B(a), B(b), and B(c) plots, and the irrigated plots were named K(a), K(b), and K(c) plots. The total area of the selected plots was 3.22 hectares. The area of the smallest plot was 0.26 hectares, and the largest one was 1.23 hectares. The total area of the rain-fed plots was 1.39 hectares with the average size of 0.46 hectares. The total area of irrigated plots was 1.83 hectares with the average size of 0.61 hectares. The average slope gradient of the study plots is 18 • ranging from 12 to 22 • . In general, terrace sizes were bigger even in the steeper slopes of rain-fed terrace as compared to less steep slope K(a). Plots K(b) and K(c) were bigger terrace size due to the lower slope gradient as compared to other plots.

UAV-Based DEM Derivation
The UAV combined with a digital camera system allows for the acquisition of usable DEM data for geomatic application such as 3D mapping and modelling [28]. In this study, the DEM data of the study area were acquired in August 2018 during the period of off farming. In this period, the irrigated fields were opened after harvesting paddy and the rain-fed cultivated lands were opened after harvesting maize and the field preparation for millet plantation begins. It avoids the disturbances of crop plants to derive high-resolution images. Photogrammetry works on the basis of two or more overlapping (stereo pair) photographs to derive the unique three-dimensional location of a set of homogeneous image points of overlapping photographs, by which the position and orientation of the camera can be computed and can be obtained a mosaic orthoimage at the end [29,30]. The first requirement is the preparation for the scaling and georeferencing for the target surface to extract oriented and scaled data. At least three or more ground control points (GCP) are required for scaling and georeferencing [31]. Altogether, 57 GCP points were fixed in six plots, and the coordinates of these GCPs were measured with total station to enable their high accuracy. GCPs were easily distinguishable from the background, which was easily marked during the processing of aerial images. The coordinates of the GCPs were used as a reference location for photogrammetric processing of the aerial images and accuracy assessment of DEM [32][33][34].
The DJI Phantom 4 Pro Plus Quadcopter Drone with Deluxe Controller-CP.PT.000549 was used mounting GO PRO camera to capture low height remote sensing imagery [35]. The flight path was defined in the UAV control software. The flight route was designed in wind direction by which flying against the wind very slowly above ground speed could be achieved. It helps to avoid motion blur because of the conventional high-flying speed. When it was flying at the low altitude, it was not possible to capture the whole area with a 60% overlap in one flight due to the low trigger frequency of the camera. The flying speed was reduced as much as possible and multiple flights at the same flight altitude were done to acquire complete dataset to compensate the gaps. The flying height was from 50 to 100 m.
Agisoft Photoscan software was used for photogrammetric process [32,36]. Feature detection algorithm was applied to align overlapping images at first. The second step was to remove inconsistencies. The third one was to apply bundle adjustment algorithm to simultaneously solve the 3D geometry of the scenes [37]. The outputs of these steps were a sparse point cloud. The GCPs on the images were identified manually, and their coordinates were imported, which were used to refine the camera calibration parameters and to optimize the geometry of the output point cloud. The GCPs coordinates were used to refine the camera calibration parameters and optimize the geometry of the output point cloud. Multi-view stereo matching algorithm was used to increase the density of sparse point cloud, which is strongly related to the number of matching key points. A polygon mesh was created from the dense point cloud and a texture map derived from all images was applied to the polygon mesh and used to create an orthophoto. Then, DEM data with pixel size of 5 × 5 cm were generated.

Soil Erosion Estimation Method
The widely used RUSLE was selected here to estimate soil erosion over the study area. This model is internationally applicable and provides comparatively reliable result [41][42][43][44]. This model is the most popular and reliable soil erosion model, which has been used around the world to predict long term rate of soil erosion from different farm size units subject to different management practices [45]. The structure of the RUSLE model can be expressed with the following equation. The average annual soil erosion rate is expressed in ton/hectare/year (t/h/y): where, E = Annual average soil erosion rate. R = Rainfall erosivity. K = Soil erodibility. LS = Topographic factors. C = Land cover management factor. P = Protection factor. According to this equation, it is known that the LS factor is the only factor that is highly connected with the surface topographic condition. For a certain plot, other factors can be thought to be invariant. Therefore, the impact from the DEM resolution can be directly assessed by calculating the LS factor with the eight DEM data. Detail processes in order to derive different factors are described.

LS Factor
The LS factor is the combined calculated values using slope gradient (s) and slope length (l). Slope gradient (s) is the steepness of slope, which has a positive impact on soil erosion rate; the slope length (l) is the distance from the point of origin of overland flow to the point where either the slope decreases to the extent that deposition begins or runoff water enters a well-defined channel. The soil erosion rate increases with the increase in slope length (l), but it has weaker strength to increase soil erosion rate with increasing slope length as compare to slope steepness (s) [45]. Topographic factors such as slope length (l), slope gradient (s), and LS factors were calculated using individual spatial resolution DEM for individual sampled plots. The following equations were used to derive slope length (l) and LS factors from DEM [46].
where s is slope gradient, and it is used in percent slope. l is the slope length in meters, which is calculated using Equation (3) A s is unit contributing area (m 2 ) [47].

R Factor
The local level climatic data are available for a very short period of time and meteorological stations are located at distant places. It is impossible to develop reliable high spatial resolution R factor data based on long term rainfall condition in the study area. Thus, rainfall erosivity (R factor) was collected from the European Soil Data Centre (ESDAC) at 1 km scale. It was calculated from the rainfall data of temporal coverage of 30-40 years with the unit of MJ mm/ha/h/y. These data were produced in 2013-2017 in collaboration between the European Commission (Joint Research Centre), the University of Basel, and the Meteorological and Environmental Institute [11]. The values of the R factor of each plot were determined according to the geolocation with the ESDAC product, with the range from 7331.9 to 9078.43 MJ mm/ha/h/y. The average value of R factor of the rain-fed plots was lower (7830.5 MJ mm/ha/h/y) than that of the irrigated plots (8528.4 MJ mm/ha/h/y). Plot K(c) has the highest value, while plot B(b) has the lowest value.

K Factor
Three soil samples from each selected plot were collected. Soil structures were identified from the field and other properties of the soil were derived by physical and chemical analysis in the laboratory. Soil erodibility factor (K) is the rate of soil loss per unit of R, as measured on a unit plot under specific condition. The soil erodibility factor (K) was calculated using the following equation: a is organic matter content in percent (%) in the soil. b is soil structure class code as defined in the U.S. Soil Survey Manual [48]. c is soil permeability class code as defined in the U.S. Soil Survey Manual [48].

C and P Factor
The land use/land cover data were derived from high resolution (5 × 5 cm) images acquired from UAV. Cultivated surfaces of the terraces were assigned as cultivated, and terrace riser and its premises were assigned as grassland. Maize is cultivated in the rain-fed outward sloping terraces in the summer season mixing with beans and vegetable crops and is fallow in the winter. The irrigated areas are paddy field with a single crop in the summer followed by wheat in the winter. Thus, the cover management factor is quite different in these cultivated plots. However, literature classifying cultivated land into the irrigated and the rain-fed are highly lacking to assign separate C and P factor. This study has borrowed the C and P factor derived from previous long-term application in Miyun reservoir area in China, because there are almost the same land use categories and crop calendar as this study [49,50]. In addition, the Miyun reservoir area has a similar seasonal rainfall pattern and altitude. Thus, the cover management factor (C) was set at 0.1 for the irrigated cultivated land, 0.35 for the rain-fed cultivated field, and 0.04 was for the grassland. The protection factor (P) was assigned 1 for all land cover categories.

Comparison Scheme
With the value setting of the RUSLE factors, the RUSLE models have been applied eight times to each plot separately, by using the topographic factors (LS) derived by eight different spatial resolutions DEM. The slope gradient, the slope length, and topographic factors (LS) were calculated using the eight different spatial resolutions DEM data, and accordingly, eight different soil erosion rates were estimated using different topographic factors (LS). Other factors (R, K, C, and P) were kept constant during the estimation. The comparison of the slope gradient (s), slope length (l), topographic factor (LS), and soil erosion of each individual plot was conducted to identify the sensitivity of these values to the spatial resolution of the DEM data.
A total percentage change was the value of total change divided by the value derived by 5 × 5 cm pixel size DEM multiplied by 100. The percentage change from one-pixel size DEM to the succeeding pixel size DEM was the percentage change in the total change.
Then, the data of individual plots were compared, then summarized into irrigated and rain-fed plots and, finally, summarized into a single data for individual changes in slope gradient, slope length, topographic factor, and soil erosion with changes in the pixel size of DEM.

Slope Gradient and Slope Length
As shown in the Figure 2a, there is a decreasing trend of slope gradient with the increase in pixel size of DEM data in rain-fed plots. The plots B(a) and B(c) exhibit similar changing pattern, but the change rate is a slightly more in plot B(c) as compared to plot B(a). Plot B(b) has the lowest change rate. A total percentage change was the value of total change divided by the value derived by 5 × 5 cm pixel size DEM multiplied by 100. The percentage change from one-pixel size DEM to the succeeding pixel size DEM was the percentage change in the total change. Then, the data of individual plots were compared, then summarized into irrigated and rain-fed plots and, finally, summarized into a single data for individual changes in slope gradient, slope length, topographic factor, and soil erosion with changes in the pixel size of DEM.

Slope Gradient and Slope Length
As shown in the Figure 2a, there is a decreasing trend of slope gradient with the increase in pixel size of DEM data in rain-fed plots. The plots B(a) and B(c) exhibit similar changing pattern, but the change rate is a slightly more in plot B(c) as compared to plot B(a). Plot B(b) has the lowest change rate. Statistically, the total change in slope gradient in plot B(c) is 45.3 (73%), which is the highest, and the lowest is in plot B(b) with the difference of 34.9 (55.5%) ( Table 1). More than 71% of the total change in the average line of the rain-fed plots happens during the changes in DEM data from 5 × 5 cm to 1 × 1 m. Thus, the change rate of the slope gradient is relatively higher from 5 × 5 cm to 1 × 1 m and with a slight increase from 1 × 1 to 10 × 10 m.
Comparatively, a similar decreasing trend was detected for the irrigated plots as shown in Figure 2b. The decreasing trend of slope gradient of plot K(a) and K(c) are nearly parallel. The values of total decrease in slope gradient in K(a) and K(c) are 34.3 and 31.8, respectively. The amount of change is lower than that of plot K(b), which has a total decrease of 46.9. The average change in slope gradient of the irrigated plots is 37.7, which is a little smaller than the average change in rain-fed plots (39.6). Similarly, the rate of change in the slope gradient is higher for the estimates using DEM data with the pixel size smaller than 1 × 1 m than that of coarse resolution such as 2 × 2, 5 × 5, and 10 × 10 m. In average, the variation in slope gradient of irrigated plots with the DEM resolution from 5 × 5 cm to 1 × 1 m occupies 77.5% of the total change. Statistically, the total change in slope gradient in plot B(c) is 45.3 (73%), which is the highest, and the lowest is in plot B(b) with the difference of 34.9 (55.5%) ( Table 1). More than 71% of the total change in the average line of the rain-fed plots happens during the changes in DEM data from 5 × 5 cm to 1 × 1 m. Thus, the change rate of the slope gradient is relatively higher from 5 × 5 cm to 1 × 1 m and with a slight increase from 1 × 1 to 10 × 10 m.
Comparatively, a similar decreasing trend was detected for the irrigated plots as shown in Figure 2b. The decreasing trend of slope gradient of plot K(a) and K(c) are nearly parallel. The values of total decrease in slope gradient in K(a) and K(c) are 34.3 and 31.8, respectively. The amount of change is lower than that of plot K(b), which has a total decrease of 46.9. The average change in slope gradient of the irrigated plots is 37.7, which is a little smaller than the average change in rain-fed plots (39.6). Similarly, the rate of change in the slope gradient is higher for the estimates using DEM data with the pixel size smaller than 1 × 1 m than that of coarse resolution such as 2 × 2, 5 × 5, and 10 × 10 m. In average, the variation in slope gradient of irrigated plots with the DEM resolution from 5 × 5 cm to 1 × 1 m occupies 77.5% of the total change.
In contrast to the variation pattern of slope gradient, the changing trend of slope length shows a positive trend with the increase in DEM pixel size. Unlike the obvious differences among different plots for the slope gradient, there are no significant differences in the changing trend of slope length among rain-fed and irrigated plots (Figure 2c,d). The plot K(b) has an undulating trend of change as compared to other plots. In the irrigated plots, there is a slow increase in the average value from 5 × 5 to 20 × 20 cm, and it starts to increase rapidly until 5 × 5 m pixel size. The average increasing rate of slope length is higher between 20 × 20 cm to 5 × 5 m, which occupies 78.1% of the total change (Table 1). In addition to presenting the estimates of slope gradient and length with different DEM resolutions, a power function can be used to fit the changing pattern of slope gradient against DEM pixel size (Figure 3a). The decreasing trend is gradually stable with the increase in DEM pixel size. The slope gradients of the irrigated plots are systematically lower than those of the rain-fed plots. Comparatively, the changing trend of slope length with the increase in pixel size presents a logarithmic trend rather than a power trend (Figure 3b). The rain-fed and irrigated plots show quite a similar pattern. For both variables, the most obvious changes are observed at the pixel size, which ranges from 5 × 5 cm to 2 × 2 m.
For better understanding of this influence, an example of the slope gradient and the slope length estimates was shown in Figure 4; with the DEM data aggregated into different resolutions, the slope gradient and the slope length of the study area can be mapped. An average slope gradient was derived at 24.2 • using 5 × 5 cm pixel size of DEM, while the slope gradient decreased along with the increase in the pixel size of DEM; finally, the slope gradient was found to be only 12.2 • in the 10 × 10 m pixel size DEM. Nearly half of the slope gradient has decreased from 5 × 5 cm to 10 × 10 m pixel sizes of DEM. Thus, it is clear that there is a significant decrease on slope gradient with increase in the pixel size of DEM. The slope length is also an important factor effecting the soil erosion rate, although it has a smaller effect than the slope gradient [46]. Thus, it is also important to know how a slope length changes with changing DEM resolution. Values of derived slope length using different resolution DEM is quite inverse to the slope gradient. The slope length increases with the increase in pixel size of DEM. The average slope length was only 33 cm in 5 × 5 cm pixel size of DEM, while it measured 115 cm in 10 × 10 m pixel size, which is a nearly three and half times increase when the pixel sizes increased from 5 × 5 cm to 10 × 10 m of DEM. This is a clear indication of high sensitivity of DEM resolution to the derivation of the slope gradient, the slope length, and the topographic factor. For better understanding of this influence, an example of the slope gradient and the slope length estimates was shown in Figure 4; with the DEM data aggregated into different resolutions, the slope gradient and the slope length of the study area can be mapped. An average slope gradient was derived at 24.2° using 5 × 5 cm pixel size of DEM, while the slope gradient decreased along with the increase in the pixel size of DEM; finally, the slope gradient was found to be only 12.2° in the 10 × 10 m pixel size DEM. Nearly half of the slope gradient has decreased from 5 × 5 cm to 10 × 10 m pixel sizes of DEM. Thus, it is clear that there is a significant decrease on slope gradient with increase in the pixel size of DEM. The slope length is also an important factor effecting the soil erosion rate, although it has a smaller effect than the slope gradient [46]. Thus, it is also important to know how a slope length changes with changing DEM resolution. Values of derived slope length using different resolution DEM is quite inverse to the slope gradient. The slope length increases with the increase in pixel size of DEM. The average slope length was only 33 cm in 5 × 5 cm pixel size of DEM, while it measured 115 cm in 10 × 10 m pixel size, which is a nearly three and half times increase when the pixel sizes increased from 5 × 5 cm to 10 × 10 m of DEM. This is a clear indication of high sensitivity of DEM resolution to the derivation of the slope gradient, the slope length, and the topographic factor.

LS Factor and Soil Erosion
It is evident from above that the DEM pixel size has a positive relation with the slope length and inverse relation with slope gradient. According to the Equation (2), the LS factor shows a decreasing trend, as shown in Figure 5a,b. There are smooth decreasing trends for the estimated LS factors in all rain-fed plots. However, it is highly fluctuating in the irrigated plot K(b) due to the high fluctuation of slope length. The average changing rate of LS factor of rain fed plots is higher than that of irrigated plots.

LS Factor and Soil Erosion
It is evident from above that the DEM pixel size has a positive relation with the slope length and inverse relation with slope gradient. According to the Equation (2), the LS factor shows a decreasing trend, as shown in Figure 5a,b. There are smooth decreasing trends for the estimated LS factors in all rain-fed plots. However, it is highly fluctuating in the irrigated plot K(b) due to the high fluctuation of slope length. The average changing rate of LS factor of rain fed plots is higher than that of irrigated plots.

LS Factor and Soil Erosion
It is evident from above that the DEM pixel size has a positive relation with the slope length and inverse relation with slope gradient. According to the Equation (2), the LS factor shows a decreasing trend, as shown in Figure 5a,b. There are smooth decreasing trends for the estimated LS factors in all rain-fed plots. However, it is highly fluctuating in the irrigated plot K(b) due to the high fluctuation of slope length. The average changing rate of LS factor of rain fed plots is higher than that of irrigated plots.  The average LS factor of the rain-fed plots has decreased 88.3% when using the DEM data with spatial resolution from 5 × 5 cm to 10 × 10 m. For the irrigated plots, the change is measured at −83.4% (Table 2). There are higher shrinking rates of LS factor from 5 × 5 cm to 1 × 1 m pixel size of DEM in all plots (including rain-fed and irrigated). There is more than an 80% reduction in the value of LS factor in these pixel sizes of DEM. The estimates of the LS factor show high sensitivity to the use of high-resolution (≤1 × 1 m) DEM data than those with coarse resolution (≥1 × 1 m).
With the estimated LS factor, combined with other factors, the soil erosion of each plot was estimated with the average value at different spatial resolutions. As shown in Figure 5c,d, there is a visible decreasing trend in soil erosion with the increase in the DEM pixel size. For the rain-fed and irrigated plots, the variation of soil erosion at different DEM pixel sizes is mainly controlled by the LS factor because of similar soil properties and rainfall conditions. Statistically, more than 90% variation of the amount of soil erosion was determined by LS factor in all individual plots at 1% significance level. Therefore, the soil erosion presents a decreasing trend for all plots with the increase in DEM pixel size.
It is evident from the Table 2 that the average amount of soil erosion of the rain-fed plots is 189.7 t/h/y when using DEM data at 5 × 5 cm pixel size. The value decreases to 57.4 t/h/y at 10 × 10 m pixel size with the reduction of 132.2 t/h/y (69.72 % decreases). Comparatively, the average soil erosion of the irrigated plots is much lower as compared to the rain-fed plots due to the differences in the topographic and C factor. The average soil erosion rate of the irrigated plot is 56 t/h/y when using 5 × 5 cm DEM data, and it changes to 11 t/h/y at 10 × 10 m pixel size. The total decrease is 44.1 t/h/y, and in percentage terms, it accounts for 79.3%. The soil erosion estimates of plot K(b) reflect the impact from the fluctuation of the LS factor. For the irrigated plots, an 83% decrease was contributed by 5 × 5 cm to 1 × 1 m pixel size of DEM in total decrease, and different plots have different proportions of contribution ranging from 80 to 87%. A decreasing rate of soil erosion with the increase in DEM pixel size is higher in the irrigated plots compared to the rain-fed plots. The decreasing rate of soil erosion with the increasing DEM pixel sizes is quite higher in irrigated plots as compared to the rain-fed plots because of the higher modification of parent slopes in irrigated plots as compared to the rain-fed one.
The trend lines for the average values of the LS factor and the soil erosion against the DEM pixel size present a similar changing pattern in both the rain-fed and irrigated plots (Figure 6a,b). A power function can be used to express the relationship. Although the function values are different for the same pixel size, all functions indicate that the decrease in the values of LS factor and soil erosion mainly happens in the range of pixel size from 5 × 5 cm to 2 × 2 m with the increase in DEM pixel size. When the pixel size is larger than 1 × 1 m, changes in LS factor and soil erosion become very slow. cm to 1 × 1 m pixel size of DEM in total decrease, and different plots have different proportions of contribution ranging from 80 to 87%. A decreasing rate of soil erosion with the increase in DEM pixel size is higher in the irrigated plots compared to the rain-fed plots. The decreasing rate of soil erosion with the increasing DEM pixel sizes is quite higher in irrigated plots as compared to the rain-fed plots because of the higher modification of parent slopes in irrigated plots as compared to the rain-fed one. The trend lines for the average values of the LS factor and the soil erosion against the DEM pixel size present a similar changing pattern in both the rain-fed and irrigated plots (Figure 6a,b). A power function can be used to express the relationship. Although the function values are different for the same pixel size, all functions indicate that the decrease in the values of LS factor and soil erosion mainly happens in the range of pixel size from 5 × 5 cm to 2 × 2 m with the increase in DEM pixel size. When the pixel size is larger than 1 × 1 m, changes in LS factor and soil erosion become very slow.

Impact of DEM Resolution on Topographic Factors
The DEM is the primary data source used for depicting surface information, which provides fundamental information for wide applications in land surface processes and environmental analysis [51,52]. In this study, the impact from the spatial resolution of DEM data has clearly reflected by the estimation of the LS factor. The spatial resolution of DEM data determines the accuracy of topographic reality, which cannot be detected with coarse resolution DEM data by missing the complexities of landform and often skipping important details of land surface [53,54]. Therefore, there were significant changes in the values of the topographic factors with the change in the spatial resolution of DEM data and the computation accuracy of the LS factor, which decreased with the increase in DEM pixel size [55]. A higher accuracy of hydrological modeling was derived by using higher resolution DEM data as compared to the use of coarse resolution DEM [19,[21][22][23]. The calculation of the LS factor using the DEM data with the pixel size ranging from 1 × 1 m

Impact of DEM Resolution on Topographic Factors
The DEM is the primary data source used for depicting surface information, which provides fundamental information for wide applications in land surface processes and environmental analysis [51,52]. In this study, the impact from the spatial resolution of DEM data has clearly reflected by the estimation of the LS factor. The spatial resolution of DEM data determines the accuracy of topographic reality, which cannot be detected with coarse resolution DEM data by missing the complexities of landform and often skipping important details of land surface [53,54]. Therefore, there were significant changes in the values of the topographic factors with the change in the spatial resolution of DEM data and the computation accuracy of the LS factor, which decreased with the increase in DEM pixel size [55]. A higher accuracy of hydrological modeling was derived by using higher resolution DEM data as compared to the use of coarse resolution DEM [19,[21][22][23]. The calculation of the LS factor using the DEM data with the pixel size ranging from 1 × 1 m to 90 × 90 m over a Burnt National Park in New South Wales, Australia, had revealed that a higher resolution DEM gives closer values to the measured values of the plots. Nevertheless, freely available low-resolution DEM is not able detect micro-complexities of topographic surface such as hill slope agriculture terraces in Nepal.
The recent advent of UAV-derived images combined with computer-based powerful photogrammetric system has improved data accuracies and provides large-scale database for small area studies [28], by which it became possible to detect more accurate and detail topographic surface information of hill slope agriculture terraces of Nepal. Although the past studies used coarser resolution DEM data as compared to present study, the changing trend of slope gradient, slope length, and topographic factors were similar with the changes in DEM resolution [56,57]. However, there were no findings revealed in this study by using the DEM resolution as high as 5 × 5 cm. Meanwhile, an important phenomenon was found that the change in the rate of LS factor is quite high, while using DEM data with pixel size below 50 × 50 cm. The result indicates that this high-resolution DEM can detect micro-topographic roughness, which was impossible in the past studies.
This study presents variations over several topographic complexities of terraced agriculture plots, including both rain-fed and irrigated plots. It partly increases the reliability of the finding through the comparative study. Generally, the average slope gradient of rain-fed plots is quite high as compared to irrigated plots, because they are made at steeper slope as compared to the irrigated plots in the study area. A total change in slope gradient is higher in rain-fed plots, but inversely, the rate of change is higher in irrigated plots because of irrigated plots strictly follow the contour when making terraces, but it is roughly followed in rain-fed plots. A sharp decrease in LS factor is recovered more in the rain-fed terraces as compared to that of the irrigated plots when using DEM data with a pixel size of 5 × 5 cm and 10 × 10 cm used. It is attributed to the big amount of roughness of the rain-fed cultivated surface due to recent tillage of cultivated surface of rainfed terraces and comparatively smooth surface in irrigated terraces. The higher rate of change in topographic factors in higher resolution DEM is resulted in due to detection of micro-topographic roughness of surfaces by high resolution DEM. In general, the identification of power function shows the changing trend of slope gradient and LS factor and logarithmic changing trend of slope length against spatial resolution DEM. It is the most important findings of the present study that were made possible by UAV-derived high-resolution DEM. More importantly, the changing trends are becoming stable when the DEM resolution is above 2 × 2 m, which partly indicates that this resolution is a threshold value for topographic calculation and estimate soil erosion of narrow cultivated terraces.

Impact of DEM Resolution on Estimated Soil Loss
It is well known that rainfall and topography are the most dominant controlling factors of soil erosion. In this study area, there are no significant differences in the values of rainfall erosivity and soil erodibility. Thus, the differences of soil erosion rates in different plots are highly controlled by a topographic factor. Soil erosion has a linear relationship with the slope gradient and the slope length because the velocity of runoff increases with the increase in slope gradient and slope length. The increase in flow depth associated with the runoff accumulation at the downslope, which increases the erosion power [58,59]. Therefore, the soil erosion estimation requires reliable DEM data, and the spatial resolution is an important aspect [57]. Although previous analysis has identified the decreasing trend of soil erosion with the decrease in spatial resolution of DEM [11,17,58,60,61], this study presents a similar study with the DEM resolution as high as 5 × 5 cm and identifies the changing pattern of soil erosion rate with the variation of DEM resolution from 5 × 5 cm to 10 × 10 m.
As shown in Figure 6, the changing rates of the estimated soil erosion are quite high in the resolution ranges from 5 × 5 cm to 1 × 1 m because of the detection of microtopographic surface of high-resolution DEM and sharp diminish of micro-topographic roughness as DEM resolution decreases. The generalization of micro-topographic surface with the decrease in DEM resolution greatly reduces the estimated rates of soil erosion. The higher rate of soil erosion in rain-fed plots as compared to irrigated plots is due to the higher slope gradient. Meanwhile, the decreasing rate of soil erosion in the rain-fed terraces is bigger than those of the irrigated ones between the pixel sizes of DEM from 5 × 5 to 10 × 10 cm is due to the higher change rate of LS factor. The detection of micro level variation of soil erosion rates was made possible by high resolution DEM in narrow agriculture terraces in the hill slope. With the increase in the DEM pixel size, the DEM data cannot depict micro level variations.
In addition to the general variation pattern, the identification of the power function to express the changing trend of soil erosion against spatial resolution of DEM is the most important finding of the present study. This study revealed that the DEM data are not able to detect the micro-topographic modification of natural slopes and accurately estimate soil erosion rate when the pixel sizes of DEM are bigger than two meters. The higher accuracy of soil erosion can be derived by using the smaller pixel size of DEM but actual soil erosion process in such a micro-topographic rough surface is still unknown because there is a lack of empirical studies using such a high-resolution DEM up to 5 × 5 cm pixel size comparing with measured data of soil loss. Thus, it needs further verification with measured values.

Implication from the Estimation of Soil Erosion Rate at a High-Resolution Level
This study is not concerned to actual soil erosion rate rather than a sensitivity analysis of estimated soil erosion rates with DEM resolution change, but it is equally important to note how far the estimated soil erosion rates are near to reality. Thus, the discussion of soil erosion rates in different parts of the world in general and in Nepal in particular will be supportive information for estimation validity. Variation of soil erosion rates in different places are due the variation of different controlling factors. Among these factors rainfall and topography are the most dominant physical factors [62]. Mountainous cropland areas with high rainfall intensity in developing countries have the higher soil erosion rates as compared to developed countries due to poor land management practices [11,[63][64][65]. The study revealed that soil erosion rates are higher in traditional agriculture areas of developing countries in the world [4] than their counter parts in developed world.
In Nepal, field measurement of soil erosion rates in Likhu Khola watershed in the central part of Middle Hill region of Nepal was 2.7 to 12.9 t/h/y, where the annual rainfall was 2141 mm [13]. Similarly, the soil erosion rate was 0.6 to 14 t/h/y in outward sloping terraces in the northeast part of Kathmandu (Nagarkot) and up to 35 t/h/y soil erosions in inward sloping red soil terraces in Pokhara and Jikhu Khola with up to 20 t/h/y [66]. There are strong variations of soil erosion rate in the hill slope agricultural area due to the variation of rainfall intensity, topography, management practices, and soil types. The government data indicate that soil erosion rates in the degraded forest and agriculture land in the central part of Middle Hill region were 31.5 to 140 t/h/y [67]. Compared with the above measurements in the hill slope terraces of Nepal, the estimates of soil erosion rate in this study are reasonably higher because of the higher rainfall intensity in the study area and detection of actual land surface of narrow agriculture terraces by high resolution DEM.
In addition to measured values, there are several estimates of soil erosion rate of different parts of Nepal using medium (30 × 30 m) to coarse (90 × 90 m) resolutions DEM [41,42,44,[68][69][70]. The use of medium to coarse resolution DEM data and relatively small rainfall of the study areas induced the low amount of soil erosion estimates. In comparison, this study detected higher soil erosion rates because of the high-resolution DEM data, which can detect the micro-variability of terraces and structures of surfaces. Further-more, it revealed how sensitivity is the estimation of soil erosion to DEM resolution. The uncertainty behind previous estimates should be considered in future soil erosion assessment. Further investigation using various resolution DEM in various surface roughness areas is required for the identification of suitable DEM size for a particular area. Finally, according to the characteristic of RUSLE model, it covers only sheet and rill erosion, which do not cover the gully erosion. Protection channels are made by farmers to protect their farm plots by gully erosion during torrential rainfall. Thus, it is not so serious for the gully erosion in the study plots. However, sometimes they fail to fully protect their farmland due to various limitations. Thus, it requires further investigation including gully erosion using high resolution DEM with the verification of actual measured data.

Conclusions
The spatial resolution of topographic data is one of the most important determinants for accurately estimating soil erosion rate in the hill slope agriculture terraces by controlling the accuracy of slope gradient and slope length. Therefore, the spatial resolution of DEM data is an important issue for modeling and estimation. Currently, most of the freely available DEMs are medium (30 × 30 m) to coarse (90 × 90 m) resolutions, which are unable to detect micro-topographic features modified by farmers in the hill slope for crop cultivation.
A higher relationship between spatial resolution of DEM and topographic factors of agriculture terraces in the hill slopes of Nepal was found. A higher rate of change was found in higher resolution DEM rather than lower resolution. It indicates that higher resolution DEM is suitable to derive a real topographic factor. There were differences in the changing values of topographic factors in rain-fed terraces and irrigated ones, the trend of topographic values is found to be similar. There is a negative power trend relationship of slope gradient with spatial resolution of DEM, but slope length has positive logarithmic relationship with DEM resolution. The LS factor and the soil erosion rates also have a negative power trend relationship with the spatial resolution of DEM. The changing trend of soil erosion was highly controlled by topographic factors in the study area rather than other factors such as rainfall erosivity and C factor. Hence, there was uniform relation of the changing trend of estimated soil erosion with DEM resolution as with LS factor.
This study concluded that DEM resolution is highly sensitive in land surface analysis for geomorphic processes and hydrological modeling. The DEM resolution less than 2 × 2 m pixel size cannot detect the narrow agriculture terraces of hill slope areas, and it cannot detect the actual rate of soil erosion. In contrast, the DEM data with the pixel size from 5 × 5 to 50 × 50 cm can depict topographic factors and soil erosion rates with high variability. Behind the changing mechanism of soil erosion estimates with different DEM resolutions, this case study well indicates the importance of the DEM resolution for soil erosion estimation for hill slope agriculture terraces in Nepal. Estimated soil erosion rates of the study can be further validated with field measured data, which are the main remaining gap of this study.