Accuracy Comparison on Culvert-Modiﬁed Digital Elevation Models of DSMA and BA Methods Using ALS Point Clouds

: High-resolution digital elevation models (HR-DEMs) originating from airborne laser scanning (ALS) point clouds must be transformed into Culvert-modiﬁed DEMs for hydrological and geomorphological analysis. To produce a culvert-modiﬁed DEM, information on the locations of drainage structures (DSs) an km in Benchmark DS data a standard assess the performance of the DSMA method compared to the BA method. A consistent methodological framework is formulated to obtain a culvert-modiﬁed DEM using DS data, using the DSMA and resultant culvert-modiﬁed DEM is then compared with BA method respectively. The DSs found from the culvert-modiﬁed DEMs were reported as true positive (TP), false positive (FP), and false negative (FN). Based on TP, FP, and FN originating from the culvert-modiﬁed DEMs of both methods, the evaluation metrics of the false positive rate (FPR) (i.e., the commission error) and false negative rate (FNR) (i.e., the omission error) were computed. Our evaluation showed that the newly developed DS data an FPR of 0.05 federal (FHWA) and 0.12 with non-FHWA roads. The FNR with FHWA roads was 0.07, and with non-FHWA roads, 0.38. The BA method an FPR of 0.28 FHWA roads and 0.62 non-FHWA the FNR BA roads based on the FPR and FNR DEM was more accurate compared with the BA method, and the formulated framework for producing culvert-modiﬁed DEMs using DSMA-based DS data was robust.


Introduction
High-resolution (HR) digital elevation models (DEMs) deciphering the fundamental earth surface processes are in frequent use and obtained from three-dimensional (3D) clouds of laser point measurements, comprising dense sampling rates (e.g., 230 laser return from a 1 m 2 area) [1,2]. Airborne laser scanning (ALS) technology has shown higher horizontal and vertical accuracies in the range of 10-15 cm, deciphering the hidden topography under the forests by penetrating through canopies [3,4]. ALS point clouds-ALS ground points (e.g., bare earth laser returns) in particular-allow the creation of HR-DEMs at a userdefined ground sampling distance (GSD) (e.g., 1 m), enabling earth surface representation at several scales through several HR-DEM variants [5,6]. An HR-DEM could adopt certain variants where each variant serves a different purpose in a particular application setting [7].
For example, airborne laser scanning (ALS) ground points are interpolated to create an HR-DEM, deciphering the Earth's surface for cartographic purposes (e.g., urban planning and making contour maps) [8]. At the same time, in the application setting of hydrology and geomorphology, HR-DEM processing must be done by incorporating some additional yet important data such as drainage structures (DSs) (e.g., bridges and culverts) so that watercourses (e.g., rivers and streams) and watershed areas can be delineated accurately [9]. Such a type of HR-DEM variant is known as a hydro-enforced DEM or hydrologically correct DEM (HC-DEM) in many parts of the world [10,11] and, in particular, a culvertmodified DEM for the U.S. Geological Survey (USGS) [12]. In the rest of the paper, the term culvert-modified DEM is adopted from the USGS [12].
A common practice among ALS data vendors, including the producers of HR-DEMs, is to classify the bridges as a separate class from ALS ground points as per the Light Detection and Ranging (LiDAR) base specifications developed by the American Society of Photogrammetry and Remote Sensing (ASPRS) [13]. One of the challenges faced while classifying ALS point clouds is the lack of information on the locations of some bridges and, particularly culverts that are difficult to discern from ALS point clouds through visual inspection [14]. Yet, an ASPRS LiDAR-based classification scheme provides a separate class code representing culverts in the classified ALS data. Subsequently, the classified ALS point clouds lack a separate culvert classification in the final deliverable of classified ALS point clouds. The locations of DSs, particularly culverts, are frequently unknown to ALS data vendors. To incorporate the information on the locations of bridges and culverts for a certain locality, the DSs are either mapped or acquired from other sources (e.g., departments of transportation (DOTs), satellite, aerial imagery, and global navigational satellite system (GNSS)-based field surveys) [10,12,[15][16][17]. Substantial information on the locations of DSs is available from DOTs because DSs fall under the jurisdiction of DOTs. Bridges and culverts are an integral part of transport infrastructures (e.g., roads and railways), enabling the unobstructed passageway of watercourses (e.g., a river or a stream). Therefore, the regular inspection, functional assessment, and repair of DSs are fundamental practices for sustainable transportation and protecting the surrounding habitat from fragmentation caused by blocked or failed DSs [18][19][20][21]. By default, the DS dataset is deemed available from DOTs. Nonetheless, the locations of a large number of DSs are unknown to transport agencies, particularly culverts [22][23][24]. For that reason, the lack of DS datasets leads to efforts in establishing digital DS databases in a state-wide capacity [25,26]. To enrich the digital databases with locations of DSs, preponderance resources are being reported (e.g., satellite and aerial imagery), as are geographical information system (GIS) techniques, and GNSS-based field surveys [27][28][29].
Yet, such databases are available for certain counties or states in the U.S. with public access. Additionally, the state-wide DS digital databases established by DOTs are considered incomplete because several roads and streets are located in private establishments, and access to those roads is out of the jurisdiction of DOTs, making it challenging to obtain extensive DS data [24]. From a global perspective, such databases are rare [10,16].

State of the Present Research in Earth Science (ES)
The increased use of ALS point clouds leads to the development of new methods comprising automated workflows to extract topographic features (e.g., roads, bridges, and buildings) [30,31]. In particular, the U.S. Geological Survey (USGS) and others have developed automated workflows to address hydrology and geomorphology in the widespread adoption of automated geomorphic feature extraction (e.g., a vector flow network of surface watercourses using a culvert-modified DEM) [12,32,33] (see Figure 1). The availability of DSs (e.g., the input culvert file) is the most essential data for obtaining a culvert-modified DEM and afterward for extracting geomorphic features correctly (see Figure 1). In their research to assess the automated workflow as shown in Figure 1, the DSs were mapped through a combination of techniques comprised of manual digitizing from aerial imagery and GNSS-based field surveys [12]. Nevertheless, ALS data consumers, particularly those involved in R&D in the earth science (ES) disciplines of hydrology and geomorphology, acquire the information on the locations of DSs from external sources [12,15,16].

Algorithmic Solutions
For some regions, external DS data sources [12,15,16] are scarce, and for that reason, algorithmic solutions were developed [10,17,35]. Among those methods, the first attempt was made to burn mapped watercourses in an HR-DEM through a process known as stream burning to create a culvert-modified DEM [36]. The stream burning approach lowers the elevation of elevated road surfaces for those locations where watercourses cross the roads. On the contrary, this approach has proven to be limited because watercourse maps are difficult to acquire.
The most frequent and standard procedure being used today is the breaching algorithm (BA) method [37] initially developed by the USGS (2010) [35] and recently further improved by Lindsay et al. (2015) [38] (For details, see Section 3.2).
However, the authors have documented the performance of the BA method to create a culvert-modified DEM. Their analysis only documented that the bridges, culverts, and roadside ditches were successfully resolved using a 1 m HR-DEM originating from ALS data. They reported an overall accuracy of about 80% [38] of DSs resolved correctly from an HR-DEM using the BA method. On the contrary, the number of DSs incorrectly resolved by their method were not documented (i.e., the commission error rate). Additionally, the number of untreated DSs (i.e., the omission error) is not explained. Therefore, the BA method requires further evaluation.

New Developments
Whether in connection with DS mapping for DOTs (see Section 1.1.1) or R & D in hydrology and geomorphology (see Section 1.1.2), the state of the present research on the subject of DS mapping reports a large variety of methods and techniques developed [12,[15][16][17]39]. Yet, no unified and automated method exists, and the need for such a unified method is reported [10,16].
Only until recently was a cost-effective and automated drainage structure mapping algorithm (DSMA) developed by Wang and Fareed (2021) [14]. The DSMA was developed using classified ALS point clouds and road centerline data. The evaluation of DS mapping was made using a benchmark DS dataset acquired from the Vermont Department of Transportation (VTrans), which was further augmented using Google Earth Street View (GE-SV) images [14]. In their research, the evaluation metrics of the prediction accuracies (e.g., precision (P), recall (R), and F1-Score) and positional accuracies, in terms of Euclidean distances of the mapped DSs compared to the benchmark DSs, were assessed. Using the prediction and positional accuracies, the authors addressed all aspects of the evaluation metrics for two different sites in Vermont, USA [14]. The assessment on DS mapping compared to the benchmark DS dataset showed that the DSMA performance was in close agreement with the benchmark DS dataset of two different landscapes (e.g., urban and rural) test sites. Therefore, the DSMA is a capable tool for addressing the DS mapping challenges in a wide area capacity. The detailed evaluation metrics of the DSMA can be found in [14].

Contributions
This paper aims at using the BA method [38] to resolve DSs from an HR-DEM to create a culvert-modified DEM (see Section 3.2). We aimed at calculating the false negative rate (FNR) (i.e., the omission error) and false positive rate (FPR) (i.e., the commission error) using the culvert-modified DEM of the BA method so that the overall performance of the BA method could be assessed, which was lacking in the past studies (see Section 1.1.3).
The DSMA is a newly developed DS mapping algorithm that needs to be quantified compared to existing standard procedures (e.g., the BA method) for obtaining a culvertmodified DEM so that the information on this topic can be updated for others (e.g., the USGS) ( Figure 1). We aimed at formulating a framework to obtain a culvert-modified DEM using the mapped DSs of the DSMA method (see Section 3.3). Afterward, using the FNR and FPR originating from the newly developed culvert-modified DEM of the DSMA method [14], the accuracy of the proposed solution will be assessed.
Finally, the FNR and FPR originating from the culvert-modified DEM of the BA and DSMA methods will be compared to assess the performance of both methods.

Study Area and ALS Survey
The focus of the present study was an area of 36 km 2 located at the latitude and longitude of 72 • 48 7.84" W and 44 • 54 41.76" N in the Missiquoi Watershed in northern Vermont, USA ( Figure 2). The land-use land-cover (LULC) matrices derived from the LULC map obtained from the Vermont Center of Geographic Information (VCGI) [40] shows that the dominant land covers were agriculture with 48% and forest with 47.5% of the 36 km 2 . The remaining part comprised of bare earth, impervious surface area (ISA), and water bodies, which represented 4.5% of the total area. The airborne laser scanning (ALS) survey in the focus study area was performed by the US Geological Survey (USGS) in the spring, beginning on 15 April 2009 and ending on 16 April 2009, with no snow cover on the ground using a Leica ALS50II laser scanner [41]. The Cessna 206s fixed-wing aircraft equipped with an ALS50II laser scanner was flown at an altitude of 1372 m a.s.l. with a nominal side overlap of 36%, comprising 158 flight lines. The Leica ALS50II is capable of collecting four laser returns per outgoing pulse at a maximum frequency of 150 kHz with a beam divergence of 0.22 mrad. The acquired ALS data was further processed into the ground and nonground points with a point density of 2 pt./m 2 . The 26 benchmark points over different topographies (e.g., bare earth, grass, and ISA) were used for ALS accuracy assessment, which showed an Root Mean Square Error (RMSE) of 8.2 cm. Afterward, the classified ALS point clouds were made available at OpenTopography by the USGS for public use [42,43]. An HR-DEM obtained from the acquired ALS data deciphering the topography of the focus area is shown in Figure 2a.

Benchmark DS Dataset Along with FHWA and Non-FHWA Roads
The cartographic representation of the road centers in vector line format, known as road centerlines data, was outsourced from the VCGI public portal [40]. The road centerlines acquired from the VCGI were administrated by federal highway authorities (FHWA), herein referred to as FHWA roads (see black lines in Figure 2b). The residential and private roads out of FHWA jurisdiction were acquired using the ISA [44] (see yellow lines in Figure 2a).
FHWA authorities collected the DS dataset (e.g., bridges and culverts (red dots)) along with FHWA roads (black lines in Figure 2a) [45]. The survey was based on GNSS coordinates of bridges and culverts, as shown with red dots in Figure 2a. In the focus study area, FHWA authorities mapped 148 DSs, largely comprised of culverts and a few bridges. The GNSS survey points of the DSs were converted into a straight line representation, as shown with the red line in Figure 2b. The straight line representation of a DS carried the information of the DS locations and their orientations under roads as shown in Figure 2b.
For the non-FHWA roads (yellow lines in Figure 2a) and those out of the jurisdiction of the state highway authorities (FHWA), a total of 18 DSs were located and digitized using Google Earth Street View (GE-SV) high-resolution panoramic images [46]. A detailed description of the acquired benchmark DS dataset could also be accessed [14].

Methods
The workflow of the proposed research is presented in Figure 3. The proposed research was carried out in four stages. (a) The ALS data was preprocessed to obtain an HR-DEM. Furthermore, FHWA and non-FHWA roads were combined into single-road centerline data at the preprocessing stage. (b) Then, the BA method was applied to the HR-DEM for creating a culvert-modified DEM. (c) The DSMA method was applied to the HR-DEM for creating a culvert-modified DEM. (d) The total number of DSs in terms of true positive (TP), false positive (FP), and false negative (FN), originating from the culvert-modified DEMs of (b) and (c), was calculated using the benchmark DS dataset as a standard reference (see Section 2.2). Finally, the TP, FP, FN values were used to obtain the evaluation metrics of the FNR and FPR originating from both methods. The details about each processing stage are provided in the following subsections.

Preprocessing
To obtain a culvert-modified DEM using the BA and DSMA methods, an HR-DEM was first developed from the ALS ground points at a GSD of 0.5 m (see Figure 2a) using the natural neighbor interpolation (NNI) technique (see Figure 3a). The choice of NNI was made because NNI was implemented in the development of the DSMA [14]. Therefore, the interpolation methods needed to be consistent for the BA [38] and DSMA [14] methods so that the results from both methods cloud be compared with certainty.
Following the HR-DEM creation, FHWA and non-FHWA roads were merged into single-road centerline data. The FHWA roads were the roads obtained from government sources, while non-FHWA roads originated from the ISA [40].

Breaching Algorithm (BA) Method
A consistent narrative reported that sinks alongside roads in an HR-DEM were usually caused by unprocessed DSs (e.g., bridges and culverts) not supplied at the HR-DEM processing stage. The researcher found that by draining those sinks with trenching conduits through elevated roads, the impact of unprocessed DSs left in the HR-DEM could be minimized [10,17]. To date, sink processing (e.g., with the breaching algorithm (BA)) is used as a standard procedure, which has been developed and improved by many researchers [35,38,47] to produce culvert-modified DEMs when comprehensive external DS data is nonexistent.
A sink is a cell or group of DEM cells with the lowest elevation compared with all immediate neighboring cells, thus making a DEM discontinuous surface at the locations of DSs (black cells in Figure 4) [48]. The breaching algorithm (BA) method was designed to process those sinks of an HR-DEM [38]. The BA method first begins by locating all sinks in a given DEM, as shown with a schematic diagram in Figure 4. Once all the sinks are found in the DEM (e.g., the black cells in Figure 4), the BA method selects one sink from the list at a time (i.e., a source sink that needs to be processed) (see the red cells in Figure 4). This is followed by finding the candidate sink(s) in the neighborhood within a user-defined search radius (red circle in Figure 4). Among the candidate sinks (green cells in Figure 4), the closest candidate is selected based on the Euclidean distances (d1, d2, and d3) computed from the center of the source sink (red) to the centers of the candidate sinks (green), as indicated by the black lines in Figure 4. Finally, the elevations of the DEM normal cell(s) within the source sink and corresponding candidate sink are lowered in such a way that the steepest descent path can be computed from a sink cell at a higher elevation to a sink cell at a lower elevation, as indicated with the yellow cells in Figure 4. The BA method iteratively fixes most of the sinks in an HR-DEM to create a culvert-modified DEM [38]. Based on the principle of connecting two sinks in close immediacy by a user-defined search radius (Figure 4), the BA method was successful in constructing line conduits through the elevated road surfaces, connecting a sink at one side of the road to a sink on another side of a road (Figure 5a). Afterward, watercourses can pass freely through newly developed line conduits of the culvert-modified DEM (Figure 5b) [38].
Using the preprocessed HR-DEM (Figure 3a), the BA method was applied to the HR-DEM to obtain a culvert-modified DEM using open-source geographical information system (GIS) software called Whitebox Geospatial Analysis Tools (GAT) [49]. The maximum user-defined search radius of 100 m was used for finding a source sink and corresponding candidate sink (see the schematic in Figure 4). The choice of a 100 m search radius ensured that all of the DSs in the study area (Figure 2a) of lengths of 100 m or less could be resolved by the BA method. This met our requirement, as the maximum length of a DS found in our study area was 70 m. To empower the BA method to resolve all the DSs from the HR-DEM to create a culvert-modified DEM (Figure 2a), an additional 30 m margin was allocated for the user-defined search radius for those DSs, in case they were overlooked by authors and larger than 70 m. This is one of the limitations of the BA method, as the user must define a user-defined search radius. For example, a DS shown in Figure 5c is about 40 m long, and it was successfully resolved by setting the user-defined search radius to 100 m. For a lower user-defined search radius, such as 39 m, the 40 m DS (Figure 5c) could remain unsolved. Therefore, the user-defined search radius was found by trial and error in their method [38].
After creating the culvert-modified DEM using the BA method. As shown in Figure 5a, watercourses can pass freely through the newly developed trenches in the culvert-modified DEM (Figure 5b). Afterward, a drainage network (DN) (i.e., watercourses map) was extracted (Figure 5b) from the culvert-modified DEM of the BA method. To find the total DSs resolved in the culvert-modified DEM. The DN was finally intersected with road centerlines (red dots in Figure 5c).

DS Mapping Algorithm (DSMA) Method
Wang and Fareed (2021) presented a new and automated approach to map DSs by developing a drainage structure mapping algorithm (DSMA) using ALS point clouds and road centerline data [14]. According to their approach, first, all types of roads (e.g., government-administrated federal highway administration (FHWA)) and those out of the jurisdiction of the FHWA (i.e., non-FHWA) were eliminated from the ALS ground points using a combined road mask created from buffers of various sizes, according to the corresponding road widths of each road class (Figure 6a). First, their evaluation showed that, by removing roads from the ALS point clouds, most of the watercourses, which otherwise were obscured by elevated roads [50], could be found (Figure 6b) [14]. Following the road removal from the ALS point clouds, the remaining ALS ground points were interpolated to create a new type of DEM known as an ALS-modified DEM (ALS-mDEM) (Figure 6), which was used to derive unobstructed DNs using the D8-method and Strahler stream orders (Figure 6b).
Second, to retrieve the locations of the DSs automatically, the derived DN from the ALS-mDEM was then intersected with the road centerlines of the FHWA and non-FHWA roads (blue line in Figure 6c) [14].
The DSMA mapping performance, compared to a trusted benchmark DS dataset, showed that a large number of DSs were mapped using a DSMA comprised of culverts [14]. The DSMA final output comprised of two datasets. First, the mapped DSs in point format were already evaluated by Wang and Fareed (2021) [14]. The second mapped DN is shown in Figure 6b. Nevertheless, The DSMA s original contribution was to provide the information on the locations and orientations (see Figure 6c,d) of DSs in a wide area capacity, and the culvert-modified DEM by default was not produced by the DSMA. However, the DS data produced by the DSMA could be useful to obtain a culvert-modified DEM using mapped DSs as the input culvert file to the add culvert tool of the Geomorphology toolbox developed by the USGS (see Figure 1).
To use the mapped DS data, further processing is required because the add culvert tool requires line DS inputs. To achieve that, using the road mask [14], the DSs were obtained as line data (blue line, Figure 6c). However, the DSs in the line conduits mapped using the DSMA (Figure 6c) were not perfectly straight lines (see Figure 6c). Nevertheless, the line conduits were generalized up to 5 m to obtain a straight line conduit representing the DSs (Figure 6d) using the ArcGIS 10.8.1 generalization tool [51]. Generalization is an automated process to remove extra nodes from a vector line by making the resulting line smoother. We tried the generalization values of 1-5 m and found that a 5 m generalization value resulted in better straight line conduits (Figure 6d). Figure 6d shows that the DSs in the straight line conduits provided the required input culvert file. To consider the generalized line as a valid input compared to the benchmark DS dataset, the mapped DS inlet and DS outlet must be aligned (see blue line in Figure 6d). Finally, the DSs in the straight line conduits were used to obtain the culvert-modified DEM using the add culvert tool of the USGS developed in ArcGIS desktop 10.8.1 (Figure 3c) [12], as shown in Figure 6e.
To further extend our approach, the culvert-modified DEM was processed twice using the mapped DS of the DSMA method. First, the HR-DEM was converted into a point LAS file, and afterward, the point LAS file was used to map DSs. Using the processing framework as shown in Figure 6, a culvert-modified DEM was produced (Figure 3c). The assessment was made to see whether a DSMA was capable of mapping DSs using an HR-DEM when ALS point clouds were not available. Second, the DSMA was used to map DSs basef on original ALS point clouds. Then, the mapped DSs were used to create a culvert-modified DEM (see the red box in Figure 3c). This was the validation step to see if the DSMA performance was similar to those DSs mapped using an HR-DEM compared to original ALS point clouds (red box in Figure 3c). Finally, DNs were then extracted from the culvert-modified DEMs and intersected with road centerlines (red dots in Figure 6f) to find the total DSs found in the culvert-modified DEMs of the DSMA method (Figure 3c).

Evaluation Metrics
To construct the evaluation metrics of the FPR using Equation (1) and the FNR using Equation (2), the total DSs found from the culvert-modified DEMs of the BA and DSMA methods were classified as TP, FN, or FN (Figure 3d).
TP represents when the watercourses correctly allowed through DSs in the culvertmodified DEM of the BA method ( Figure 5c) and DSMA method (Figure 6c) compared to a corresponding benchmark DS dataset. FP represents when watercourses were allowed through artificially constructed DSs compared to the benchmark DS datasets. Finally, FN represents when watercourses were not allowed through DSs compared to the benchmark DS dataset:

DS Mapping Assessment of the BA Method Using the HR-DEM
We first report the results of the descriptive statistics of DSs as TP, FP, and FN, obtained from a culvert-modified DEM with the BA method (see Table 1). The DS classification was carried out separately, along with FHWA and non-FHWA roads (see Section 2.2), to assess the performance rate of the BA method along with each road type. We are the first who have investigated the causes and effects of FN and FP cases compared to the benchmark DS dataset for the BA method in order to assess its performance. Therefore, the FN and FP cases were further categorized into three distinct cases based on our evaluation. Case-1 and Case-2 were associated with FNs, and Case-3 was associated with FPs. The detailed descriptive statistics are shown in Table 1, whereas each case scenario is discussed in Section 4.3, compared to the DSMA method using benchmark DS data as a standard reference.

DS Mapping Assessment of the DSMA Method Using an HR-DEM and ALS Point Clouds
The number of DSs found from culvert-modified DEM using a DSMA-based DS dataset are also reported as TP, FP, and FN in Table 2 for comparison with the BA method (Table 1).  Table 2 shows that the total number of DSs classified as TP, FP, and FN from the culvert-modified DEM based on the mapped DSs using a HR-DEM equaled the TP, FP, and FN from the culvert-modified DEM based on mapped DSs using ALS point clouds only (see Table 2). Table 2 shows the interesting trends which emerged, as the culvert-modified DEM obtained either HR-DEM-based mapped DS data or the original ALS point clouds were not affected by the data source (see Table 2).
The following evaluation emphasizes that, alone, the HR-DEM is equally useful compared to ALS point clouds for mapping DSs using the DSMA method and thereafter producing the culvert-modified DEM. Therefore, the DS could be mapped by transforming the HR-DEM into point data without compromising the DS mapping results.
The accuracy evaluation when creating the culvert-modified DEM using DS data originating from the DSMA method (Table 2) compared to the BA method (Table 1) showed interesting trends. Table 2, compared to Table 1, shows that the FPR of 0.05 observed for the FHWA roads and 0.12 observed for the non-FHWA roads were much lower compared with the BA method, which was 0.28 for the FHWA roads and 0.62 for the non-FHWA roads. Similarly, the FNR of 0.07 for the FHWA roads and 0.38 for the non-FHWA roads were found, compared with the BA method, which found 0.32 for the FHWA roads and 0.61 for the non-FHWA roads (for comparisons, see Tables 1 and 2).
In the context of creating the culvert-modified DEM, the emerging trend of higher evaluation matrices of the DSMA method compared to the BA method was further investigated and compared based on each case scenario of Table 1.

Comparative Assessment of BA and DSMA Methods
Our evaluation suggests that the DSMA is based on a better approach than the BA method, resulting in a more precise culvert-modified DEM. The DSMA first corrects the ALS point clouds from obstructing man-made features (e.g., roads) before it can generate an ALS-modified DEM (ALS-mDEM), and afterward, a DN, as shown in Figure 6c [14]. By removing the man-made impediments from the ALS-mDEM, the stage is set to map almost all watercourse pathways (Figure 6c) passing through roads. Bridges and culverts are mostly found over watercourses, and by intersecting with road centerlines, the found DSs are approximated very close to the benchmark DS dataset, which makes the DSMA approach more robust to map DSs. Afterward, the mapped DSs could be incorporated in an HR-DEM using the formulated framework (see Figures 3c and 6).
On the contrary, the BA method is for solving the sinks by trenching the conduits based on the Euclidean distance between two sinks very close to each other (see Figure 4). The culvert-modified DEM was obtained using the BA method, but not all of the DSs were processed. Our investigation showed that by following the underlying principle of fixing sinks caused by DSs, the BA method resulted in three different case scenarios (Table 1) where it lacked resolution for several DSs in order to create a culvert-modified DEM.

Case-1: No Source or Candidate Sink Found
To trench a conduit through elevated roads, the BA method must find an upstream sink (i.e., a source sink) and a corresponding downstream sink (i.e., the candidate sink) (see Figure 4). Figure 7a shows a DS deciphered correctly (TP) in a culvert-modified DEM (blue line) and the same DSs deciphered by the DSMA method (yellow line) compared to the benchmark DS (red line).
However, this is not always the case for several DSs where an upstream source sink was found while a corresponding downstream candidate sink was missing. Therefore, trenching a conduit by the BA method was not possible, causing an FN (Figure 7b). This finding also validates the Duke et al. (2003) observations [50] as a DN becomes discontinuous (Figure 7b) at both sides of the road. For several other locations, the source sink and corresponding candidate sink were missing both sides of the road as shown in Figure 7c, causing FN cases for the BA method originating from the Case-1 scenario (Table 3).  On the contrary, Figure 7b,c shows that the DSMA method did not suffer through such problems associated with the BA method (e.g., searching for source and corresponding candidate sinks) ( Figure 4). Instead, the DSMA found the watercourse pathways using an ALS-mDEM, whereas the roads were already removed so that watercourses passed freely on both sides of the road. The DN was based on the flow direction and flow accumulation metrics [14], as indicated by the yellow line, deciphering all DSs mapped correctly (TP) by the DSMA method compared with the BA method (blue line).

Case-2: No Solution within the Found Sinks
It has been observed that a source sink and corresponding candidate sink both were found upsteam on the same roadside within the user-defined search radius (see Figure 4), as shown in Figure 8. Yet, the candidate sinks on the other side of the road were absent, causing no solution to be found by the BA method, resulting in FN cases ( Figure 8 and Table 1). Compared with the BA method, the DSMA method successfully found the underlying DS, as shown by the yellow line in Figure 8. Due to the above-stated reasons, the number of FN cases was always higher for the BA method for the FHWA and non-FHWA roads compared with the DSMA method (Tables 1 and 2). This led the BA method to an FNR of 0.32 along FHWA roads and 0.61 along non-FHWA roads, which was significantly higher than the DSMA method s FNR of 0.07 along FHWA roads and 0.38 along non-FHWA roads (Tables 1 and 2).

Case-3: Wrong Solutions within the Found Sinks
In Case-3, the BA method was successful in finding a source sink upstream and a corresponding candidate sink downstream ( Figure 9). Afterward, a conduit was trenched between the source sink and corresponding candidate sink. However, the proposed solution was wrong compared with the benchmark DS dataset, indicating no DS for that location and therefore causing FP cases ( Figure 9 and Table 1).
All of the FPs originating from the BA method belonged to Case-3 for the FHWA and non-FHWA roads (Table 1). The nature of the FPs caused by the DSMA algorithm was briefly explained by Wang and Fareed (2021), and readers should refer to [14], which states that the DSMA creates some artificial DN causing FP cases; therefore, that error propagates to the culvert-modified DEM based on mapped DSs using the DSMA method. However, the total FP cases caused by the DSMA method was much lower compared with the BA method (see Tables 1 and 2).
The overall comparison of FP and FN cases caused by the BA method showed similar statistics, where the FP cases were almost equal to the FN cases (Table 1). Therefore, the FPR and FNR for the BA method were similar along the FHWA and non-FHWA roads, respectively.
To understand this trend, the case scenarios shown in Figures 7-9 were further analyzed. We found that for the BA method, if a right solution was not found (FN) (see Figures 7 and 8), a corresponding false solution (FP) was found ( Figure 9). Therefore, if the watercourses aere not allowed through the DS (FN), the watercourses were then allowed through artificial DSs (Figure 9) created by the BA method (FP).
For that reason, the FPR and FNR were similar in the focus study area. For FHWA roads, an FPR of 0.28 was reported, corresponding to an FNR of 0.32. Similarly, for non-FHWA roads, we found an FPR of 0.62 with a corresponding FNR of 0.61. The BA method, as reported by others, can fix major bridges and culverts [40]. Therefore, the FPR and FNR along non-FHWA roads were found to be higher because non-FHWA roads are frequently comprised of culverts of smaller sizes, as observed by Wang and Fareed (2021) [14].
Compared with the BA method, the DSMA approach showed reasonably consistent performance along FHWA and non-FHWA roads [14]. For FHWA roads, an FPR of 0.05 was found, corresponding to an FNR of 0.07, which was significantly lowere than the BA method. Similarly, for non-FHWA roads, we observed an FPR of 0.12, corresponding to an FNR of 0.38.

Discussion
In the context of creating a culvert-modified DEM, this study found that the DSMA method was useful for mapping DSs in a straight line conduit format automatically, which could be incorporated in an HR-DEM to obtain a culvert-modified DEM (see Figure 1), compared with manual methods, although with available information on the locations of DSs (e.g., GNSS locations), a significant amount of time is required to digitize the DS locations into straight line conduits manually for producing culvert-modified DEMs. However, the DSMA supplies the DSs in straight line conduits automatically, reducing the manual digitization workload.
Compared with the BA method, the DSMA showed improved performance for creating culvert-modified DEMs. Our evaluation showed that DSMA-based DSs resolved a much higher number of DSs from HR-DEMs compared with the BA method (Tables 1 and 2). Figures 7-9 show that the culvert-modified DEM obtained using DS data mapped by the DSMA always allowed the watercourses to connect both sides of the roads compared with the BA method, which does not permit the watercourses to pass either side of the road [50,52].
However, for the DSMA method, the watercourses did not exactly pass through the DSs compared with the benchmark DS dataset (see Figures 7-9). The detailed statistics on the positional accuracy of the mapped DS using the DSMA are shown in Table 3 and were found by Wang and Nadeem (2021) [14]. The positional accuracy was based on the Euclidean distances between the mapped DSs and their corresponding benchmark DS reference data. Thus, positional accuracy determined how well the DSs were mapped by the DSMA compared to the standard reference (e.g., the benchmark DS dataset). A detailed description of the positional accuracy of mapped DSs was given in [14].
The descriptive statistics, as shown in Table 3, indicate that the positional accuracy was about 10 m for the median and third quartile of the mapped DSs by the DSMA. The positional accuracy in the range of 10 m was acceptable according to Lidberg et al. (2017) [39], where the positional offset of watercourses was set to 10 m in their research. Nevertheless, for exceptionally fewer cases, a maximum offset of 72.53 m was expected for those regions where the topography is nearly flat (e.g., slope ≈ 0 • ) [14]. For such regions, once DSs are mapped using the DSMA method, it is possible to adjust the mapped DS to their original positions manually using hillshade images, orthophotos, and GE-SV images [53][54][55]. Such required manual interventions are exceptionally few and only associated with flat topographies [14].
Our evaluation showed that the DSMA-based mapped DSs were useful in the following application environments.
First, we have shown that for those regions having no access to ALS point clouds but an available HR-DEM, the HR-DEM could be used to map DSs using the DSMA method. This is useful for those who have access to ALS derivatives such as an HR-DEM and access to original ALS point clouds that is constrained by many factors (e.g., copyrights and data privacy). Second, to derive a culvert-modified DEM, extensive ground surveys could besignificantly reduced, as the DSMA provides extensive DS maps in a wide area capacity using either an HR-DEM or is based on ALS point clouds.

Conclusions
A newly developed DSMA s (https://github.com/Nadeem1geo/Drainage-Structure-Mapping-Algorithm-DSMA-, accessed on 6 February 2021) performance was assessed in the context of producing a culvert-modified DEM compared to the breaching algorithm (BA) method. The culvert-modified DEM is essential in the ES disciplines of hydrology and geomorphology. Previously, the breaching algorithm was commonly used to create culvert-modified DEMs. While the DSMA is a newly developed algorithm that has never been tested in this context, both methods were assessed using a benchmark DS dataset along FHWA and non-FHWA roads. The classified ALS point clouds and road centerlines were acquired for the study area of 36 km 2 in Vermont, USA.
The evaluation of prediction accuracies of both methods was carried out using the FPR and FNR with the benchmark DS dataset. We found that the DSMA was robust and created a more precise culvert-modified DEM compared with the BA method by returning the lowest FPR of 0.05 along FHWA roads and 0.12 along non-FHWA roads. Similarly, the FNR of 0.07 for the FHWA roads and 0.38 along the non-FHWA roads were observed for the DSMA method, making prediction accuracies closer to the benchmark DS dataset. In contrast to the DSMA, the BA method showed an FPR of 0.28 along the FHWA roads and 0.62 along the non-FHWA roads. Similarly, the BA method resulted in an FNR of 0.32 along the FHWA roads and 0.61 along the non-FHWA roads, making prediction accuracies five times lower than the DSMA method. We conclude that the DSMA can potentially be used to develop a culvert-modified DEM, compared with the BA method, and the methodological framework is presented in this study.
Author Contributions: Conceptualization, Nadeem Fareed and Chi-Kuei Wang; methodology, Nadeem Fareed; software, Nadeem Fareed; formal analysis, Nadeem Fareed; writing-original draft preparation, Nadeem Fareed; writing-review and editing, Chi-Kuei Wang and Nadeem Fareed; visualization, Nadeem Fareed; supervision, Chi-Kuei Wang. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement:
The data that support the findings of this study can be requested from the corresponding author.