Detection of Catchment-Scale Gully-Affected Areas Using Unmanned Aerial Vehicle ( UAV ) on the Chinese Loess Plateau

The Chinese Loess Plateau suffers from serious gully erosion induced by natural and human causes. Gully-affected areas detection is the basic work in this region for gully erosion assessment and monitoring. For the first time, an unmanned aerial vehicle (UAV) was applied to extract gully features in this region. Two typical catchments in Changwu and Ansai were selected to represent loess tableland and loess hilly regions, respectively. A high-powered quadrocopter (md4-1000) equipped with a non-metric camera was used for image acquisition. InPho and MapMatrix were applied for semi-automatic workflow including aerial triangulation and model generation. Based on the stereo-imaging and the ground control points, the highly detailed digital elevation models (DEMs) and ortho-mosaics were generated. Subsequently, an object-based approach combined with the random forest classifier was designed to detect gully-affected areas. Two experiments were conducted to investigate the influences of segmentation strategy and feature selection. Results showed that vertical and horizontal root-mean-square errors were below 0.5 and 0.2 m, respectively, which were ideal for the Loess Plateau region. The overall extraction accuracy in Changwu and Ansai achieved was 84.62% and 86.46%, respectively, which indicated the potential of the proposed workflow for extracting gully features. This study demonstrated that UAV can bridge the gap between field measurement and satellite-based remote sensing, obtaining a balance in resolution and efficiency for catchment-scale gully erosion research.


Introduction
Soil erosion is a serious environmental problem which causes high economic costs [1].Many approaches have been proposed for monitoring and predicting soil erosion at different scales [2][3][4][5].Gully erosion, the main type of soil erosion by water, is a major land degradation process that adversely affects land management and agriculture [6].In the past century, an increasing number of researches focus on the gully erosion due to its ubiquity and severity [7].Evaluation of gully-affected areas is the basis for control and monitoring of gully erosion; as such, detection of these areas has become a growing interest in gully erosion community [8][9][10].
Gully extraction method depends on the development of geographic data acquisition.Field investigation is the most traditional method.In the early stage, tapes, rulers, and micro-topographic profilers are commonly used [11,12].Recently, latest measurement technologies, such as terrestrial laser scanning (TLS) [13], airborne laser scanning (ALS) [14], 3D photo-reconstruction [15], and total station [16], are adopted.Both TLS and ALS are able to efficiently acquire high resolution dataset for erosion monitoring and related fluvial geomorphology studies [3,17,18].Recent experiments in different regions including Australia, Peru, and England prove that the TLS gains more adoption than ALS for catchment-scale studies because it is more flexible and accurate [19][20][21].However, the heavy fieldwork complicates the application of such method for a large region [18].
With a large number of earth observation satellites launched, abundant imagery can be used to assess gully erosion [22].Optical satellite sensors, including the medium-resolution imagery (e.g., Landsat and Spot) and high-resolution imagery (e.g., IKONOS and Quick Bird), are increasingly available.The current frequently used imageries in geoscience are Worldview-3 and Pleiades [23,24], with a ground sampling distance in panchromatic mode of 0.31 and 0.5 m, and 1.00 and 2.00 m in multispectral mode, respectively.Remote sensing based method for extracting gully features shows two advantages over field investigation.On the one hand, multi-temporal and multi-resolution geographic data covering almost the world are easily obtained [25][26][27].On the other hand, remote sensing not only provides spectrum and texture information but also generates the digital elevation models (DEMs) when the satellite possesses stereoscopic capabilities [28,29].Remote sensing is mainly carried out on satellite or manned aircraft, which are good options for regional-scale data acquisition because of their wide coverage and stable performance [30].Nevertheless, the application of conventional platforms is limited for the increasing demanding of environmental modeling at catchment scale because of their high cost, low flexibility, and poor spatial and temporal resolution [31].
The use of unmanned aerial vehicle (UAV) can bridge the gap between field investigation and satellite and aircraft-based remote sensing [32].UAV-based remote sensing is suitable for catchment-scale surveys, which is more flexible than traditional remote sensing.Recently, many publications prove that UAV can be regarded as a credible tool for monitoring soil erosion [33], coastal area [34], precision agriculture [35], glacier dynamic [36], and landslide [37], in which the created high-resolution DEM and ortho-image mosaics can provide further detailed information.
As the method for data acquisition transfers from the field investigation to remote sensing, the extraction method also improves greatly.In the early stage, visual interpretation, which is based on the spectrum difference and interpreter's knowledge, is the main choice [38,39].However, manual interpretation is restricted to low efficiency and uncertainty which has been replaced by automatic method.Pixel and object-based are two kinds of automatic method for gully feature extraction.Although, the pixel-based method has been applied in many studies [40,41], the object-based method is proven to be more advanced because it can integrate the spectral, shape, and textural information instead of spectral information only [42,43].The recent papers show that object-based method is the mainstream in processing high-resolution imagery [44][45][46][47].
Although UAV shows potential in gully features extraction at catchment scale, existing studies are conducted in limited regions, such as Morocco [32] and Saxon loess province [33].According to experiments all over the world, the contribution of gully erosion to overall soil loss rates and sediment production rates by water erosion range from 10% to 94% [6].Significant differences in gully erosion conditions lead to varying gully sizes, shapes, and densities in different regions.Thus, further studies are needed to investigate UAV-based method for gully feature extraction in other typical areas.
The Chinese Loess Plateau is known for its serious soil erosion and land degradation, which are caused by both natural and anthropogenic factors [48].Gully erosion in this region accounts for 60%-70% of the total soil loess [49], with a large number of developed gullies shaping the distinct loess landform [50].The Loess shoulder-lines divide the total area into upland and gully-affected areas with totally different topographic signatures, land use, and types of soil erosion [51].Gullies in this region can be divided into three types: floor, bank and hillslope gullies [52].Although some studies have discussed the extraction of gully features in Chinese Loess Plateau, specifically on catchment scale [53][54][55], the use of UAV is currently limited.The present study aims to: (1) apply the UAV on the Chinese Loess Plateau to generate the high-resolution DEMs and ortho-mosaics; and (2) investigate the object-based method for detection of gully-affected areas by using UAV-acquired dataset.

Study Area
The Loess Plateau is located in the middle and upper reaches of the Yellow River, and covers an area of 640,000 km 2 [56].Based on the loess landscape, Loess Plateau can be divided into loess tableland and loess hilly regions [51].Loess tableland can be regarded as the early stage of loess landform, which retains its original paleolandforms with a large plain area.With the development of erosion, the plain area will be reduced gradually and replaced by the narrow loess ridges and isolated loess mounds [57].The loess hilly regions suffer from more intensive soil erosion than loess tableland region; consequently, the geomorphologic landscape of the former becomes ruined and complex.
In this case, two study sites were selected to represent two typical loess landforms.The study site representing loess tableland is located in the northwest of Changwu County, which is a part of Xialiu catchment.This study site, which consists of plain area, hillslope area and deep gullies, covers approximately 2.33 km 2 .The elevation in this site ranges from 981 to 1220 m.The other study site is a part of Zhifang catchment, covering 3.42 km 2 , and located in southern Ansai County.This site is the typical loess hilly region with the elevation ranging from 1129 to 1417 m.The gully-affected area covers more than 60% of total area, with many intensely developed gullies.The locations and pictures of study sites are shown in Figure 1.

Method
The two methods employed in this paper (Figure 2) are: (1) UAV photogrammetry for obtaining high-resolution DEM and ortho-mosaics; and (2) object-based detection of gully-affected areas by using generated datasets.During data acquisition, strict regulations should be observed from outdoor survey to indoor image processing so that the accuracy metrics can satisfy the requirements.For detection part, objected-based approach and random forest (RF) classifier were applied with two experiments for optimizing segmentation and classification steps.

Method
The two methods employed in this paper (Figure 2) are: (1) UAV photogrammetry for obtaining high-resolution DEM and ortho-mosaics; and (2) object-based detection of gully-affected areas by using generated datasets.During data acquisition, strict regulations should be observed from outdoor survey to indoor image processing so that the accuracy metrics can satisfy the requirements.For detection part, objected-based approach and random forest (RF) classifier were applied with two experiments for optimizing segmentation and classification steps.

Method
The two methods employed in this paper (Figure 2) are: (1) UAV photogrammetry for obtaining high-resolution DEM and ortho-mosaics; and (2) object-based detection of gully-affected areas by using generated datasets.During data acquisition, strict regulations should be observed from outdoor survey to indoor image processing so that the accuracy metrics can satisfy the requirements.For detection part, objected-based approach and random forest (RF) classifier were applied with two experiments for optimizing segmentation and classification steps.

UAV Description
In this study, microdrone md4-1000 is applied (Figure 3a), which is a miniaturized VTOL aircraft (Vertical Take off and Landing).The entire system consists of an aerial vehicle, a ground station with the software for mission planning and flight control, a radio control transmitter, and a telemetry system [58].Md4-1000 can fly for approximately 45 min (depending on payload and wind) at 15 m/s.This microdrone can fly by remote control or automatically using microdrone GPS waypoint navigation software.
The maximum recommended payload of md4-1000 UAV is 0.80 kg.A digital system camera, Sony ILCE-7R (Sony Corporation, Tokyo, Japan), was mounted on the UAV.This camera acquired images in true color (Red, Green, and Blue bands) with 8-bit radiometric resolution.The camera was also equipped with a 50 mm zoom lens.The camera's sensor is 7360 × 4912 pixels, and the images are stored in a secure digital SD-card.Image triggering is activated by the UAV according to the programmed flight route.The camera was also equipped with GPS receiver, altimeter and wind meter, so that the onboard computer system can record a timestamp, the GPS location, flight altitude, and vehicle principal axes (pitch, roll, and heading) at the time of each shoot.In this study, microdrone md4-1000 is applied (Figure 3a), which is a miniaturized VTOL aircraft (Vertical Take off and Landing).The entire system consists of an aerial vehicle, a ground station with the software for mission planning and flight control, a radio control transmitter, and a telemetry system [58].Md4-1000 can fly for approximately 45 min (depending on payload and wind) at 15 m/s.This microdrone can fly by remote control or automatically using microdrone GPS waypoint navigation software.
The maximum recommended payload of md4-1000 UAV is 0.80 kg.A digital system camera, Sony ILCE-7R (Sony Corporation, Tokyo, Japan), was mounted on the UAV.This camera acquired images in true color (Red, Green, and Blue bands) with 8-bit radiometric resolution.The camera was also equipped with a 50 mm zoom lens.The camera's sensor is 7360 × 4912 pixels, and the images are stored in a secure digital SD-card.Image triggering is activated by the UAV according to the programmed flight route.The camera was also equipped with GPS receiver, altimeter and wind meter, so that the onboard computer system can record a timestamp, the GPS location, flight altitude, and vehicle principal axes (pitch, roll, and heading) at the time of each shoot.

Outdoor Survey
The outdoor survey can be divided into two steps: image acquisition and ground control survey.For image acquisition, three individuals, namely, a ground station operator, a radio control pilot and a visual observer, were needed for the secured use of the UAV.Weather condition should be also considered, particularly to prevent influence of the wind and rain.In this study, the operations of the UAV were performed in the morning when wind was relatively low to ensure flight stability and image quality [36].The flight plans in the two study areas were designed using the flight planning software.The images were taken with an average overlap of 70% in flight direction and 60% overlap in flight strip.It should be noted that the spacing between flight lines in Ansai is uneven.It is because that more terraces are distributed in the right part so that less flight line is needed due to the gentle topography.
Ground control points (GCPs) are required to update the horizontal and vertical accuracies.Although the direct georeferenced system is employed in some studies [59], there are high requirements for the camera and GPS which are expensive and heavy at this stage.Two kinds of GCPs, namely, base control points (BCPs) and photo control points (PCPs), were implemented to obtain a highly precise coordinate information.In each study area, three BCPs were installed for building the GPS network covering the entire study area.Each BCP was measured for nearly 1.5 h with Topcon HIPERⅡG in static mode (Figure 3b).Compared with BCP, more PCPs are needed.The placement of PCPs should be between the flight lines and located in the overlap areas.In this study, 38 PCPs and 53 PCPs were collected for Changwu and Ansai, respectively, by using the GPS rover (Figure 3c).Horizontal coordinates were referenced to China Geodetic Coordinate System 2000, and the vertical values were referred to the National Vertical Datum 1985.Locations of the GCPs and flight lines can be seen in Figure 4.

Outdoor Survey
The outdoor survey can be divided into two steps: image acquisition and ground control survey.For image acquisition, three individuals, namely, a ground station operator, a radio control pilot and a visual observer, were needed for the secured use of the UAV.Weather condition should be also considered, particularly to prevent influence of the wind and rain.In this study, the operations of the UAV were performed in the morning when wind was relatively low to ensure flight stability and image quality [36].The flight plans in the two study areas were designed using the flight planning software.The images were taken with an average overlap of 70% in flight direction and 60% overlap in flight strip.It should be noted that the spacing between flight lines in Ansai is uneven.It is because that more terraces are distributed in the right part so that less flight line is needed due to the gentle topography.
Ground control points (GCPs) are required to update the horizontal and vertical accuracies.Although the direct georeferenced system is employed in some studies [59], there are high requirements for the camera and GPS which are expensive and heavy at this stage.Two kinds of GCPs, namely, base control points (BCPs) and photo control points (PCPs), were implemented to obtain a highly precise coordinate information.In each study area, three BCPs were installed for building the GPS network covering the entire study area.Each BCP was measured for nearly 1.5 h with Topcon HIPER II G in static mode (Figure 3b).Compared with BCP, more PCPs are needed.The placement of PCPs should be between the flight lines and located in the overlap areas.In this study, 38 PCPs and 53 PCPs were collected for Changwu and Ansai, respectively, by using the GPS rover (Figure 3c).Horizontal coordinates were referenced to China Geodetic Coordinate System 2000, and the vertical values were referred to the National Vertical Datum 1985.Locations of the GCPs and flight lines can be seen in Figure 4.

Indoor Image Processing
Indoor image processing includes aerial triangulation, DEM generation, and ortho-mosaics calculation.Aerial triangulation refers to the process for obtaining the true positions and orientations of the images by using the photographs covering the same object, exposed from different positions.The number of GCPs is usually limited; as such, a large number of tie points are created for conjugate points identified across multiple images [60].The input data for aerial triangulation contain scanned images, camera calibration, and the GCPs.Aerial triangulation was performed in the Trimble's Inpho 6.0 by using the bundle-block adjustment, which is widely adopted in previous studies [31].To make triangulation successful, root mean square error (RMSE) was calculated for the orientation of image block.If the error exceeds the threshold value, refinement is needed by removing the blunder tie points.
After the triangulation process, the DEMs and ortho-mosaics can be created with a high resolution because of the stereoscopic overlap of the images [32].Only image pairs with well-suited overlap and good quality were selected.The removal of redundant images can reduce the processing time.MapMatrix 4.0 developed by Visiontek Inc. was used for this step.Although MapMatrix can produce a DEM automatically using the feature-matching technology, manual editing is also necessary for the study areas with large relief.In this study, automatic algorithm was also applied for generating a robust surface.Subsequently, the terrain features, such as the peck points, ridge lines, and valley lines, were collected for creating a high-resolution DEM.To improve the DEM quality, the

Indoor Image Processing
Indoor image processing includes aerial triangulation, DEM generation, and ortho-mosaics calculation.Aerial triangulation refers to the process for obtaining the true positions and orientations of the images by using the photographs covering the same object, exposed from different positions.The number of GCPs is usually limited; as such, a large number of tie points are created for conjugate points identified across multiple images [60].The input data for aerial triangulation contain scanned images, camera calibration, and the GCPs.Aerial triangulation was performed in the Trimble's Inpho 6.0 by using the bundle-block adjustment, which is widely adopted in previous studies [31].To make triangulation successful, root mean square error (RMSE) was calculated for the orientation of image block.If the error exceeds the threshold value, refinement is needed by removing the blunder tie points.
After the triangulation process, the DEMs and ortho-mosaics can be created with a high resolution because of the stereoscopic overlap of the images [32].Only image pairs with well-suited overlap and good quality were selected.The removal of redundant images can reduce the processing time.MapMatrix 4.0 developed by Visiontek Inc. was used for this step.Although MapMatrix can produce a DEM automatically using the feature-matching technology, manual editing is also necessary for the study areas with large relief.In this study, automatic algorithm was also applied for generating a robust surface.Subsequently, the terrain features, such as the peck points, ridge lines, and valley lines, were collected for creating a high-resolution DEM.To improve the DEM quality, the editing work is needed to modify the elevation to eliminate the influences of the buildings and vegetation.Finally, the DEM at 1 m resolution was exported, and the ortho-image mosaics at 0.2 m were calculated after input images were orthorectified.

Accuracy Assessment
To investigate the accuracy of the DEM data, some independent ground truth points are needed.The spatial continuous data are the most suitable options [61]; however, no reference data are available for the selected study area with such high-resolution requirements.Instead, some independent check points (ICPs) were surveyed, with 33 ICPs in Changwu and 37 ICPs in Ansai.According to the requirements for monitoring the gullies, the majority of these points were chosen along the gully border lines (Figure 4).The measurements of ICPs were similar to those of GCPs.The vertical error was calculated by the differences between the elevation of ICPs and the DEM data.The horizontal deviation was calculated by measuring the displacements between the ICPs and the ortho-mosaics.The RMSE was used to evaluate the accuracy according to Equation (1): where Z di is the i-th measured value from DEM or ortho-mosaics, and Z ri is the corresponding reference data from the ICPs.

Data Preparation
In addition to the high-resolution DEM and ortho-mosaics, a large amount of information should be prepared for gully feature extraction.In existing papers, topographic features derived from DEM exhibit potentials in distinguishing the natural objects related to terrain factors, such as landslide [62], riverscape [63] and gully erosion [45].In the present study, shaded relief, slope, roughness, and specific catchment area (SCA) [64], which were involved into the further segmentation and classification steps, were created from the original DEM.
Reference data are the indispensable part of this method.The reference data in training area were employed to optimize the parameters and determine the rule-sets.For the test area, reference data were used for accuracy assessment.Field investigation was conducted in April 2016 to obtain the credible reference and further understanding of the shape, distribution, and development of gullies in these two study areas.The final reference data were digitized manually via integrating imagery features, expert knowledge, and the field investigation.

Segmentation
Compared with the pixel-based method, the analysis unit of object-based method transfers from one pixel to object, in which image segmentation is key step.Among all kinds of segmentation method, multi-resolution image segmentation algorithm (MRIS) [65] is the most widely used and is implemented in the eCognition software.The MRIS is a region growing algorithm, which starts from one pixel by merging adjacent pixels based on a heterogeneity criterion [66].
Image segmentation directly influences the final classification result.The selection of layers is the key issue to obtain the suitable objects using eCognition.In existing studies for gully feature extraction, only imagery is used in the segmentation step, and the terrain information is ignored.DEMs and its derived topographic layers has proven the potentials to improve the segmentation result in geomorphic mapping [67,68].In this paper, two segmentation strategies were designed to evaluate whether the topographic information could improve the segmentation accuracy: (1) only ortho-mosaics were used with three optical bands (DOM strategy); and (2) the DEM and ortho-mosaics, including three optical bands, slope, roughness, and shaded relief (hillshade), were integrated (OSRH strategy).
The heterogeneity criterion depends on the segmented layers and parameters, including scale, shape, and compactness.Scale factor controls the object size, which is the key factor in MRIS.Two methods are commonly used to obtain the optimized scale parameter: (1) local variance based method [69]; and (2) spatial autocorrelation based method [70].Since Drȃguţ [71] proposed the Estimation of Scale Parameter (ESP) using local variance method in 2010, ESP-Tool has been widely employed to determine the scale parameter because it is user-friendly.The improved version, which supports the multiple layer parameterization, was employed in the present study [72].In addition to the scale factor, shape factor determines the weight of shape heterogeneity, and color heterogeneity should also be emphasized.Within the shape criterion, the weight assignment between compactness and smoothness is determined by compactness factor.However, more attention has been paid to the selection of scale parameter, litter discussion is made for shape and compactness.
In this study, scale, shape, and compactness were all considered for determining the most suitable parameter combination.First, the visual assessment, which can reduce the time cost, was used for determining the possible range of each parameter (Table 1).Second, ESP-Tool was used for each parameter combination.Any deleted scale value within the range could be saved for further investigation.Third, the goodness of segmentation was evaluated using the following metrics: under-segmentation, over-segmentation and Euclidian distance [73,74]: where r is the reference data and s is the segmentation dataset.OR and UR are area-based metrics which describe the match degree between ground-true and segmentation results.ED is the combined metrics which can be regarded as the "closeness" to the ideal segmentation result.

Image Classification: RF
RF is an ensemble machine learning algorithm proposed by Leo Breiman [75].This supervised method can build a forest of decision tree based on the classification and regression tree algorithm (CART).Compared with CART algorithm, two improvements were achieved by using bootstrap sampling and randomized node optimization.For each tree, the bootstrap, a kind of random sampling with replacement, was performed; hence, the same sample may be selected several times, while others could not be selected at all.To estimate the RF model performances, the out-of-bag error was calculated in the cross-validation technique.For each node, the decision was selected with user-defined number of feature selected randomly, instead of using all the features.Those two mechanisms caused RF to exhibit high classification ability and processing speed.Currently, RF is widely used in earth science community [76].
Two user-defined factors are involved in RF model: the number of trees that will grow (nTree) and the number of randomly selected features (mTry).Based on the results published to data, nTree is typically set to 500, and mTry is set to the square root of the variable number [76].Currently, many studies focus on the influence of feature selection on the classification accuracy.In RF method, variable importance (VI) is employed to optimize the feature space based on the Mean Decrease in Gini (MDG) or the Mean Decrease in Accuracy (MDA).MDG describes how each variable contributes to the homogeneity of nodes.MDA is a measure based on the out-of-bag error calculation.The more accuracy decreases because of the exclusion of a variable, the higher MDA it will be.In this study, MDA was used for feature selection as the majority reports did.
Metric calculation of each object is needed between the segmentation and classification steps.Four kinds of metrics, including spectral, textural, geometric, and topographic information, were calculated per object (Table 2).Among the total 36 variables, spectral and geometric features are commonly used in object-based image analysis.Texture features, which can improve the classification, are also widely used.Eight features derived from Grey Level Co-occurrence Matrix (GLCM) proposed by Haralick [77] were calculated based on two bands.One feature is the red band, which represents the image texture, and the other is shaded relief layer created in ArcMap which represents the terrain texture.In addition, topography is a key factor for the initiation and development of gullies [78].In this paper, the mean values of the slope, roughness, shaded relief, and SCA were calculated for each object.

DEM and Ortho-Mosaics Generation
After the field investigation and indoor image processing in April 2016, the high-resolution DEM and ortho-mosaics were generated.Figure 5 shows the shaded relief images of both areas.The 1 m resolution DEM can reflect the overall erosion conditions, morphological characteristics, and micro topography with detailed information.Two enlarged areas for each study areas were selected representing the natural gully-affected areas and man-made terrace areas.For gully erosion research, the DEM-based digital terrain analysis is important.The ortho-mosaics, in which the image context and 3D perception are useful for visual interpretation and landscape extraction, are also needed.As shown in Figure 6, a large amount of land surface information, such as the cropland, vegetation, and houses, can be expressed.
representing the natural gully-affected areas and man-made terrace areas.For gully erosion research, the DEM-based digital terrain analysis is important.The ortho-mosaics, in which the image context and 3D perception are useful for visual interpretation and landscape extraction, are also needed.As shown in Figure 6, a large amount of land surface information, such as the cropland, vegetation, and houses, can be expressed.Accuracy assessment is necessary before using data.The ICPs for both areas were measured during the outdoor survey.The vertical error was calculated between the heights of ICPs and corresponding values from DEM.Furthermore, the horizontal error was determined according to the differences between ICPs and the ortho-mosaics.Figure 7 and Table 3 summarize the statistics.The vertical RMSEs in Changwu and Ansai are 0.245 and 0.339 respectively.The horizontal errors are relative lower, with RMSE reaching 0.083 and 0.143.Many factors, including the GPS system error, measurement error for GCPs and the errors caused by indoor image processing, contribute to the output error.The terrain characteristics exert considerable influence on the accuracy.The two selected study areas are all located in Loess Plateau.The large topographic relief causes difficulty during GCPs surveying.In addition, the photogrammetric restitution and DEM extraction in these two study areas need further manual work, which also contributes to the error estimate.Hence, the accuracy is relatively lower than that of UAV photogrammetry in coastal areas [34] and mudflat [79].The reason why Changwu owns higher accuracies than Ansai can also be explained by this reason.representing the natural gully-affected areas and man-made terrace areas.For gully erosion research, the DEM-based digital terrain analysis is important.The ortho-mosaics, in which the image context and 3D perception are useful for visual interpretation and landscape extraction, are also needed.As shown in Figure 6, a large amount of land surface information, such as the cropland, vegetation, and houses, can be expressed.Accuracy assessment is necessary before using data.The ICPs for both areas were measured during the outdoor survey.The vertical error was calculated between the heights of ICPs and corresponding values from DEM.Furthermore, the horizontal error was determined according to the differences between ICPs and the ortho-mosaics.Figure 7 and Table 3 summarize the statistics.The vertical RMSEs in Changwu and Ansai are 0.245 and 0.339 respectively.The horizontal errors are relative lower, with RMSE reaching 0.083 and 0.143.Many factors, including the GPS system error, measurement error for GCPs and the errors caused by indoor image processing, contribute to the output error.The terrain characteristics exert considerable influence on the accuracy.The two selected study areas are all located in Loess Plateau.The large topographic relief causes difficulty during GCPs surveying.In addition, the photogrammetric restitution and DEM extraction in these two study areas need further manual work, which also contributes to the error estimate.Hence, the accuracy is relatively lower than that of UAV photogrammetry in coastal areas [34] and mudflat [79].The reason why Changwu owns higher accuracies than Ansai can also be explained by this reason.Accuracy assessment is necessary before using data.The ICPs for both areas were measured during the outdoor survey.The vertical error was calculated between the heights of ICPs and corresponding values from DEM.Furthermore, the horizontal error was determined according to the differences between ICPs and the ortho-mosaics.Figure 7 and Table 3 summarize the statistics.The vertical RMSEs in Changwu and Ansai are 0.245 and 0.339 respectively.The horizontal errors are relative lower, with RMSE reaching 0.083 and 0.143.Many factors, including the GPS system error, measurement error for GCPs and the errors caused by indoor image processing, contribute to the output error.The terrain characteristics exert considerable influence on the accuracy.The two selected study areas are all located in Loess Plateau.The large topographic relief causes difficulty during GCPs surveying.In addition, the photogrammetric restitution and DEM extraction in these two study areas need further manual work, which also contributes to the error estimate.Hence, the accuracy is relatively lower than that of UAV photogrammetry in coastal areas [34] and mudflat [79].The reason why Changwu owns higher accuracies than Ansai can also be explained by this reason.Nevertheless, in the present study, the accuracy remains satisfactory according to the requirements for detection of gully-affected areas.Nevertheless, in the present study, the accuracy remains satisfactory according to the requirements for detection of gully-affected areas.

Detection of Gully-Affected Areas
According to the discussion in Section 3.2.2,two segmentation strategies were used.The goodness assessments of segmentation using the candidate parameters are shown from Tables 4-7 to determine the most suitable segmentation parameters.ED value is used for determining the most desirable segmentation parameters.The error metrics based on OSRH strategy are obviously lower than the corresponding values based on DOM strategy, which prove that topographic information can improve the segmentation results significantly.Two sets of parameters based on OSRH were finally obtained: 481, 0.4, and 0.7 for Changwu, and 327, 0.3, and 0.7 for Ansai.

Detection of Gully-Affected Areas
According to the discussion in Section 3.2.2,two segmentation strategies were used.The goodness assessments of segmentation using the candidate parameters are shown from Tables 4-7 to determine the most suitable segmentation parameters.ED value is used for determining the most desirable segmentation parameters.The error metrics based on OSRH strategy are obviously lower than the corresponding values based on DOM strategy, which prove that topographic information can improve the segmentation results significantly.Two sets of parameters based on OSRH were finally obtained: 481, 0.4, and 0.7 for Changwu, and 327, 0.3, and 0.7 for Ansai.According to the segmentation results, RF model can be used for predicting the gully-affected area.The 15 selected features in the two study areas are shown in Figure 8 after calculating the mean value of MDA from 10-fold replicate runs.The topographic features in the two study areas rank on the top, thereby indicating that the gully distributions are closely related to the characteristics of terrain.Texture features also dominate the VI ranking, with seven selected in Changwu and eight selected in Ansai.Notably, the texture features derived from DEM achieve relatively better positions in VI rank than those of ortho-mosaics.This is because the terrain texture reflects the land surface without the influences of vegetation and man-made buildings; hence, terrain texture possesses better performance for improving the accuracy than image texture.Except brightness, which ranked third place in Ansai, almost all of the spectral and shape metrics rank beyond 15th place, thereby exhibiting minimal effect in improving the classification accuracy.To investigate the influences of feature number on the classification accuracy, the RF model was built using different number of features based on the VI rank.The classification assessments of the RF model are shown in Figures 9 and 10.In general, both Changwu and Ansai achieve their highest F-measures using all features, with 84.62% in Changwu and 86.46% in Ansai.The accuracies in Changwu appear stable when the selected features are more than 8.When the feature number becomes less than 8, the producer's accuracy (PA) decreases dramatically leading a low F-measure.However, the user's accuracy (UA) in this range even experiences a slight increase.Such behavior of RF classifier is caused by the class imbalance problem, in which the minority class will be over predicted [80].The estimation of balance class would be regarded as an intelligent guess on designing the RF model [62].In this paper, a further discussion about such problem is not conducted only because it is not serious within the recommended range .The results in Ansai are almost the same.The classification results can be acceptable if the selected number is more than 10.When the feature number is less than 10, PA and UA decrease simultaneously.The class imbalance problem does not happen in Ansai, which can be explained as the natural class-distribution in Ansai is more close to the balanced sample than Changwu.The RF models based on all the features display the most suitable performances; thus, these models were used for obtaining the final gully-affected areas in both areas (Figure 11).To investigate the influences of feature number on the classification accuracy, the RF model was built using different number of features based on the VI rank.The classification assessments of the RF model are shown in Figures 9 and 10.In general, both Changwu and Ansai achieve their highest F-measures using all features, with 84.62% in Changwu and 86.46% in Ansai.The accuracies in Changwu appear stable when the selected features are more than 8.When the feature number becomes less than 8, the producer's accuracy (PA) decreases dramatically leading a low F-measure.However, the user's accuracy (UA) in this range even experiences a slight increase.Such behavior of RF classifier is caused by the class imbalance problem, in which the minority class will be over predicted [80].The estimation of balance class would be regarded as an intelligent guess on designing the RF model [62].In this paper, a further discussion about such problem is not conducted only because it is not serious within the recommended range .The results in Ansai are almost the same.The classification results can be acceptable if the selected number is more than 10.When the feature number is less than 10, PA and UA decrease simultaneously.The class imbalance problem does not happen in Ansai, which can be explained as the natural class-distribution in Ansai is more close to the balanced sample than Changwu.The RF models based on all the features display the most suitable performances; thus, these models were used for obtaining the final gully-affected areas in both areas (Figure 11).The Ansai area suffers from more serious gully erosion with steeper topography, but it achieves a higher accuracy than Changwu.This phenomenon can be explained by the sample strategy for generating training and test dataset.An integrated watershed could be regarded as the ideal unit to generate the representative training samples.In this study, the experiments in Ansai were in accordance with this requirement.The left sub-watershed was used as the training area, and the right part, which included another sub-watershed and a sector of hillslope area, was selected as the test area.For Changwu, the whole survey area only contains one watershed; hence no natural boundary was present for generating the training dataset.The current sample strategy may be biased, which leads to a lower accuracy than expected.

Comparison with Existing Studies
UAV photogrammetry, a cutting-edge technology, was applied in two catchments of Loess Plateau of China.The fieldwork of UAV-based method is much less than the conventional surveying method; meanwhile, the operation flexibility is higher than established RS-based method.
As the demand for high-resolution topographic data increases dramatically, UAV-based photogrammetry is undergoing explosive development.Tarolli stated that high resolution topographic data can provide an opportunity for better understanding the earth surface process [81].Hence, the UAV photogrammetry possesses various application fields besides the gully features extraction present in this paper.The Ansai area suffers from more serious gully erosion with steeper topography, but it achieves a higher accuracy than Changwu.This phenomenon can be explained by the sample strategy for generating training and test dataset.An integrated watershed could be regarded as the ideal unit to generate the representative training samples.In this study, the experiments in Ansai were in accordance with this requirement.The left sub-watershed was used as the training area, and the right part, which included another sub-watershed and a sector of hillslope area, was selected as the test area.For Changwu, the whole survey area only contains one watershed; hence no natural boundary was present for generating the training dataset.The current sample strategy may be biased, which leads to a lower accuracy than expected.

Comparison with Existing Studies
UAV photogrammetry, a cutting-edge technology, was applied in two catchments of Loess Plateau of China.The fieldwork of UAV-based method is much less than the conventional surveying method; meanwhile, the operation flexibility is higher than established RS-based method.
As the demand for high-resolution topographic data increases dramatically, UAV-based photogrammetry is undergoing explosive development.Tarolli stated that high resolution topographic data can provide an opportunity for better understanding the earth surface process [81].Hence, the UAV photogrammetry possesses various application fields besides the gully features extraction present in this paper.
The object-based approach is also used in this paper, which has proven its ability in different regions.Nonetheless, only few studies adopt this method in the Loess Plateau region for extracting gully features.In previous studies, significant amount of work focuses on the extraction of gully border line (also named loess shoulder line).Although the gully border line can also be used for determining the gully-affected areas, the extraction methods are uncertain and not universal [82,83].The accessibility of data exerts significant influence on the method design.In the existing studies, high-resolution DEM is not available; hence, extracting the border line based on the gray difference is more practicable than the area-wide mapping.The proposed method exhibits several advantages.First, the information derived from the high resolution DEM and ortho-mosaics is fully used, which is in line with the development trend.Second, our method is a general workflow, which can be expanded to other catchments easily.Third, the detection accuracies are high and stable.In addition, the object-based approach can support hierarchical classification; consequently, both gully-affected areas and independent gullies can be extracted.

Limitations for the Application of UAV in Loess Hilly Region
The time cost for the use of UAV in two study areas is more than expected, about two days were spent for outdoor survey in Ansai and one and half day in Changwu.The placing and surveying of GCPs, specifically for the GCPs located in the gully-affected areas which are hard to reach, are time consuming.Although photogrammetric processing is sensible to the distribution and number of GCPs [32], a balance between the accuracy and time cost should be considered.Such balance is a common problem for the area characterized with steep topography.For indoor processing, the feature-based matching technique adopted in Inpho and MapMatrix allows the automatic generation of the DEM and ortho-mosaics.However, automatic processing without operator interventions may cause considerable systematic errors in areas with steep topography.Therefore, a semi-automatic workflow was applied in our study areas.These two reasons reduce the efficiency of UAV photogrammetry in Loess Plateau compared with the published studies in areas with gentle topography.
This study mainly aims to detect gully-affected area, which can be regarded as a specific type of landform [84].Some studies focused on extracting specific types of gullies [85].In the Chinese Loess plateau region, gullies can be divided into three kinds: hillslope, bank and floor gullies.As shown in Figures 5 and 6, the 1 m DEM is sufficient for detecting gully-affected area, and hillslope and bank gullies can also be recognized.UAV photogrammetry can be used for gully feature mapping in different levels.The existing studies showed its potentials in hillslope gully extraction [33,70].However, challenges still exist.The floor gullies, developing in the floor of valley, are difficult to be reflected using UAV photogrammetry because of the deep location and small terrain openness; the traditional field measurement is still more suitable in this case.For bank gullies, the UAV photogrammetry possesses difficulties for modeling the gullies in considerably steep slope, due to the limitation of shooting angle.The aerial oblique imaging which supplements the traditional vertical photography can be used in the further study to overcome such restriction [72,86].

Conclusions
Inspired by the application of UAV in environmental modeling, we applied a low-cost and highly efficient device for gully feature extraction in the Chinese Loess Plateau region.Although the topographic relief in the two selected study areas is high caused by the gully erosion, the UAV can still obtain the 1 m DEM and 0.2 m ortho-image mosaics successfully.According to the highly detailed productions, an object-based method combined with the RF classifier was used to detect the gully-affected areas.
Three contributions were provided by this paper: (1) an integral workflow for detection of gully-affected areas that includes UAV photogrammetry for data generation and object-based approach for feature mapping; (2) the segmentation step was optimized by considering the improvements of topographic features; and (3) the experiments in the two catchments of Chinese Loess Plateau provided a new method for extracting gully features in this region.
Further study is required to increase the accuracy for detecting gully-affected areas and independent gullies.In addition, multi-temporal UAV data will be acquired to monitor gully development.

Figure 1 .
Figure 1.Location of study areas, the worldview-3 images and the photos captured during the filed investigation in (A) Ansai and (B) Changwu

Figure 1 .
Figure 1.Location of study areas, the worldview-3 images and the photos captured during the filed investigation in (A) Ansai and (B) Changwu.

Figure 2 .
Figure 2. Flowchart of the proposed methodology.

Figure 2 .
Figure 2. Flowchart of the proposed methodology.

Figure 3 .
Figure 3. Illustration of the unmanned aerial vehicle (UAV) and outdoor survey: (a) microdrone md4-1000; (b) base control points surveying; and (c) the rover used to measure photo control points.

Figure 3 .
Figure 3. Illustration of the unmanned aerial vehicle (UAV) and outdoor survey: (a) microdrone md4-1000; (b) base control points surveying; and (c) the rover used to measure photo control points.

Figure 4 .
Figure 4. Overview of the survey design including the flight lines, locations of base control points (BCPs), photos control points (PCPs) and independent check points (ICPs).The background images are acquired from Google TM Earth.Transformations from China Geodetic Coordinate System 2000 to Mercator were made for all ground control points (GCPs) and fight lines using the projection tool in ArcMap.

Figure 4 .
Figure 4. Overview of the survey design including the flight lines, locations of base control points (BCPs), photos control points (PCPs) and independent check points (ICPs).The background images are acquired from Google TM Earth.Transformations from China Geodetic Coordinate System 2000 to Mercator were made for all ground control points (GCPs) and fight lines using the projection tool in ArcMap.

Figure 7 .
Figure 7. Boxplots of error measured between independent control points and the generated DEM (Vertical) and ortho-mosaics (Horizontal).

Figure 7 .
Figure 7. Boxplots of error measured between independent control points and the generated DEM (Vertical) and ortho-mosaics (Horizontal).

Figure 9 .
Figure 9. Evolution of classification accuracy related to the number of selected variables in Changwu.

Figure 10 .
Figure 10.Evolution of classification accuracy related to the number of selected variables in Ansai.

Figure 9 . 20 Figure 9 .
Figure 9. Evolution of classification accuracy related to the number of selected variables in Changwu.

Figure 10 .
Figure 10.Evolution of classification accuracy related to the number of selected variables in Ansai.Figure 10.Evolution of classification accuracy related to the number of selected variables in Ansai.

Figure 10 .
Figure 10.Evolution of classification accuracy related to the number of selected variables in Ansai.Figure 10.Evolution of classification accuracy related to the number of selected variables in Ansai.

Table 1 .
The selection range of parameter for segmentation.

Table 2 .
Overview of the features adopted for gully-affected areas detection in this paper.

Table 3 .
Statistics of the vertical error and horizontal error (all values in m).

Table 3 .
Statistics of the vertical error and horizontal error (all values in m).

Table 5 .
Segmentation accuracy metrics for Changwu based on OSRH strategy.

Table 6 .
Segmentation accuracy metrics for Ansai based on DOM strategy.

Table 7 .
Segmentation accuracy metrics for Ansai based on OSRH strategy.