Anadem: A Digital Terrain Model for South America

1. Introduction

Digital representations of the Earth’s surface provide essential data for a wide range of applications, enhancing our understanding of hydrology, landscapes and terrains, and are key components in many studies, such as topographic mapping [1,2], flood modeling [3,4,5], landslide analysis [6,7], and atmospheric boundary theory [8].

Digital elevation models (DEMs) can be generated on large scales from different sources, including cartographic data, optical satellite imagery, radar satellite data, and airborne data. At present, satellite data is most commonly used for acquiring digital representations of the Earth’s surface, primarily due to its cost-effectiveness and accuracy [9]. As a result, numerous DEMs based on optical and radar satellite data have been made publicly available, with fully or near-global coverage.

The Shuttle Radar Topography Mission (SRTM) was a pioneering initiative providing the first high-resolution DEM with almost global coverage [10]. SRTM was acquired using a C-radar instrument (31 mm wavelength), which was responsible for populating global topographic information [11]. Subsequently, satellite missions have continued to make significant improvements in terms of the accuracy and spatial resolution of DEMs. Examples include the Advanced Spaceborne Thermal Emission and Reflection Radiometer Global DEM (ASTER GDEM) [12], the Advanced Land Observing Satellite (ALOS) World 3D 30 meters (AW3D30) [13], and Copernicus GLO 30 (COPDEM30) [14].

The accuracy of DEMs is significantly influenced by the presence of vegetation [8,15,16], primarily because the most commonly used remote sensing acquisition techniques, such as interferometry (e.g., radar sensors, using X-radar instruments) and stereometry (e.g., panchromatic sensors), are unable to capture elevation information for the bare earth surface [17,18]. As a result, numerous studies have sought to develop solutions to mitigate the impact of vegetation bias on global DEMs and to create digital terrain models (DTMs) [19,20,21,22,23]. Baugh et al. [19] addressed this issue by effectively removing 50–60% of the estimated vegetation height in the SRTM, thereby improving its hydrodynamic modeling accuracy. Applying a similar approach, O’Loughlin et al. [21] incorporated remote sensing products, including vegetation density and height data, into different empirical equations to determine the most effective means of minimizing vegetation bias in the SRTM dataset.

Light detection and ranging data (LiDAR) are not affected by vegetation. Therefore, the altitude measured by sensors such as the Global Ecosystem Dynamics Investigation (GEDI) and ICESat (Ice, Cloud, and Land Elevation Satellite) can be considered as the altitude of the terrain. However, because it is not an imaging sensor, information from LiDAR does not have global coverage. For this reason, LiDAR data can only be used to directly correct DEMs in a few pixels. To correct the effect of altitude on most pixels, the solution that has been adopted is to develop an empirical model that relates the height of the vegetation, or more precisely, the magnitude of the error that the vegetation causes in DEMs [21,22,23].

For example, Yamazaki et al. [22] introduced The Multi-Error-Removed Improved-Terrain (MERIT) model, which utilized global vegetation databases and LiDar data from ICESat to reduce vegetation bias in the SRTM model. Similarly, Magruder et al. [23] utilized the differences between SRTM elevations and ICESat-2 terrain heights, canopy cover, and terrain slope to reproduce and remove vegetation biases from the SRTM model.

As the availability of remote sensing data increases and further development of empirical models takes place, the associated increase in complexity may require new approaches. In this context, machine learning advancements have significantly enhanced the utilization of extensive data to create new topographic datasets based on existing remote sensing sources [20,24,25,26,27]. Consequently, machine learning has emerged as a powerful tool for generating high-accuracy DTMs comprising large amounts of remote sensing data, including altimetry measurements acquired using LiDAR instruments, such as ICESat and GEDI. For example, Liu et al. [28] used the Random Forest (RF) model to enhance the accuracy of MERIT, resulting in an 18% improvement in flood simulation relative to the original dataset. Additionally, Hawker et al. [20] introduced the Forest and Buildings Removed Copernicus DEM (FABDEM), which mitigates biases from vegetation and urban buildings by applying a random forest model and incorporating altimetry data from the GEDI dataset.

Despite some advancements in methodologies to remove vegetation bias in global DEMs, most studies still rely on SRTM, despite the fact that it is outperformed by more recent DEMs [29]. Moreover, the potential of machine learning methods to remove vegetation bias is yet to be fully investigated; many such methods are not yet openly available. This study aims to propose a new DTM for South America, named ANADEM – Agência Nacional de Águas DEM (acronym of National Water Agency DEM), by removing bias vegetation from COPDEM. First, we developed a methodology to remove vegetation influence using machine learning and remote sensing datasets. Subsequently, we validated our DTM by comparison with terrain elevation measurements acquired using the ICESat-2 mission and available global DEMs.

ANADEM was developed for applications by the Brazilian water agency and offers several advantages, notably its calibration, which is specifically designed to improve its accuracy over South America. Additionally, we expect that, through using open-source initiatives and freely available data, we will enhance its use in applications requiring topographic information, ensuring easy accessibility to DTM data for South America.

2. Methods and Data

2.1. Methods Overview

The model workflow in this study was developed to facilitate applicability to different study areas, requiring only terrain elevation and multispectral data to correct vegetation bias from the selected DEM (COPDEM). We applied the model to produce a DTM open data for South America (ANADEM), using machine learning and multispectral information from Landsat and Sentinel, and elevation profile from GEDI. An overview of the methodology is shown in Figure 1, with further details provided in the subsequent sections.

2.2. Study Area

ANADEM was developed for the South America continent, which covers an area of 17.8 million km² and is characterized by diverse climate and land cover conditions (Figure 2). The region is primarily covered by tropical forests (47%), extensive areas of pasture, croplands, grasslands, and savannas. The tree structure in South America exhibits substantial variability in canopy density (sparse to dense) and heights, which range from short (0 to 2 meters) to tall (15 to 35 meters).

2.3. COPDEM

The Copernicus GLO-30 (COPDEM) [14] is a resampled version (with a 30-meter spatial resolution) derived from the WorldDEM dataset, which is based on radar interferometry data (X-band) acquired by the TanDEM-X mission [31] between 2011 and 2015. The creation of the WorldDEM dataset also involved filling data gaps, incorporating data from other global DEMs, and applying filtering processes. Additional filtering was performed over water bodies by employing a flattening process that ensured a more consistent hydrological network. Recent studies have demonstrated that COPDEM is the most precise freely available global DEM [29,32,33,34], motivating its use in this study.

2.4. GEDI

The Global Ecosystem Dynamics Investigation (GEDI) is a LiDAR instrument mounted on the International Space Station since 2019, with a footprint diameter of 25 meters and an acquisition data interval of 60 meters [35]. GEDI data provide valuable information on canopy structure by measuring the time difference between the emitted and return pulse, which allows the representation of vegetation structure by identifying signals from the ground and from the vegetation canopy. The main objective of the GEDI mission is to provide high-accuracy and high global coverage data on canopy structure, which can be used to estimate biomass [36,37,38] and canopy height [28,39,40].

In this study, GEDI data were acquired from NASA’s Earth Data website (https://search.earthdata.nasa.gov/). Elevations were referenced to the World Geodetic System 1984 (WGS84) ellipsoid and height metrics from the L2A product. Then, we derived the relative height (RH) metrics for each measurement, representing the cumulative energy between the lowest detectable return (ground elevation) and the highest detectable return (canopy top). To ensure data quality, we applied a filtering process to exclude lower quality GEDI data points.

2.5. Multispectral Data

Multispectral data available between 2013 and 2020 from Landsat and Sentinel were used to compute vegetation indices, which were then integrated into the classifier to capture the vegetation characteristics of the landscape (see Section 3.2). To ensure data quality, a filtering process was implemented to select only images with less than 20% cloud cover.

The Landsat 8 satellite (launched in 2013) was used due to its compatibility with the period of elevation data acquisition from COPDEM. We also selected the highest quality data globally available, from collection 2, which incorporates atmospherically corrected land surface reflectance. The satellite provides high-resolution information (30–100 meters) of the Earth’s surface at 16-day intervals, captured by the Operational Land Imager (OLI) and Thermal Infrared Sensor (TIRS) instruments.

Additionally, data from the Sentinel-2 mission, acquired using a pair of satellites equipped with the MultiSpectral Instrument (MSI), were selected. The spatial resolution of the Sentinel bands ranges from 10 to 60 meters. Both satellites operate with a five-day interval between acquisitions.

2.6. ICESat-2

ICESat-2 (Ice, Cloud and land Elevation Satellite) was launched in 2018, marking the second generation of satellite mission with the purpose of measuring ice, atmospheric and topographic parameters via LiDAR measurements [41]. ICESAT-2 is equipped with the Advanced Topographic Laser Altimeter System (ATLAS), which has a footprint diameter of 17 meters, and a spacing of 90 meters between measurements. We obtained data between January 2020 and December 2020, applying the ICESat-2 ATL03 product (https://icesat-2.gsfc.nasa.gov/) to validate our methodology based on independent data. ICESat-2 ATL03 incorporates essential information, such as quality flags and atmospheric and georeferenced corrections, and is also referenced to the WGS84 ellipsoid.

2.7. Data Composition

Landsat and Sentinel multispectral data were used to compute spectral indices, which were subsequently used as inputs for our model. Spectral indices aggregate spatial information related to vegetation, which generally improves the final results of the classifier. After conducting preliminary analyses of different spectral indices, we selected three specific indices, namely the Normalized Difference Moisture Index (NDMI), the Modified Soil Adjusted Vegetation Index (MSAVI), and the Enhanced Vegetation Index (EVI), which demonstrated superior discrimination of vegetation and effectively aggregated information into the model (Supplementary Material S1). The detailed equations for each spectral index are provided in Table 1.

2.8. Vegetation Bias Removal Algorithm

To remove the vegetation bias from the COPDEM product, we used a non-parametric machine learning algorithm, the gradient tree boost classifier [45], which is available in the Google Earth Engine (GEE) cloud platform [46]. The gradient tree boost classifier [45] is based on tree decision learning methods and has several advantages over the traditional random forest model due to its optimization techniques (boosting) aimed at minimizing errors in the prediction process. When using this method, each new tree is generated with consideration of the residual errors of the previous tree. Additionally, the gradient tree boost method demonstrated greater robustness that the traditional random forest model in generating accurate results with low sensitivity to outliers and a reduced processing time.

Initially, we divided the study area into 93 (5° × 5°) tiles to avoid processing errors and improve the geographical representativeness of vegetation characteristics. Next, we trained the classifier by randomly selecting 50,000 samples of GEDI points for each tile and extracting the satellite indices values for each point. We defined the machine learning classifier parameters by testing different optimization settings and selected those with the optimal applicability, highest accuracy, and shortest processing times (Supplementary Material S1). We then applied the gradient tree boost classifier with 200 trees and selected the Huber statistics as a loss function.

Results from the classifier were post-processed to correct exaggerated vegetation removal and minimize noise. We utilized the fractional vegetation cover (fvc) mosaic generated with Sentinel-2 multispectral data to ensure that locations with sparse vegetation underwent lower vegetation removal than dense vegetation. Water body elevation values remained equal to those used by COPDEM, and an additional step was applied for locations adjacent to water bodies in order to avoid lowered values and ensure hydrological consistency. After subtracting the vegetation bias estimated from the COPDEM, a smoothing focal median filter with a 3 × 3 pixel was applied to reduce noise.

2.9. Performance Assessment

To assess the accuracy of the proposed DTM using independent terrain elevation data, we used elevation terrain information obtained by the ICESat-2 satellite. Information from ICESat-2 is available in the ‘hf_best_fit’ band from the ATL03 product, which provides the best fit ground elevation estimative. We generated a random sample of 20,000 points from ICESat-2 for each tile, which were then compared against other global digital elevation and terrain models, including SRTM, COPDEM, MERIT, and FABDEM.

The following metrics were calculated to enable a statistical assessment of each model evaluated: the root mean squared deviation (RMSE) between the predicted and observed values (Equation 1), and the mean vias error (BIAS), which is the sum of the differences between the observed and predicted values (Equation 2).

$$\mathit{R}\mathit{M}\mathit{S}\mathit{E}=\surd \left(\frac{1}{\mathit{n}}\right)\sum _{\mathit{i}}^{\mathit{n}}{(\mathit{o}\mathit{b}\mathit{s}\mathit{e}\mathit{r}\mathit{v}\mathit{e}\mathit{d}-\mathit{m}\mathit{o}\mathit{d}\mathit{e}\mathit{l})}^2$$

(1)

$$\mathit{B}\mathit{I}\mathit{A}\mathit{S}=\frac{1}{\mathit{n}}\sum _{\mathit{i}}^{\mathit{n}}(\mathit{o}\mathit{b}\mathit{s}\mathit{e}\mathit{r}\mathit{v}\mathit{e}\mathit{d}-\mathit{m}\mathit{o}\mathit{d}\mathit{e}\mathit{l})$$

(2)

3. Results and Discussion

3.1. Model Validation

The validation results for vegetation bias removal in the COPDEM model are presented as a new DTM model named ANADEM. Overall, ANADEM led to an 85.6% reduction in vegetation bias relative to COPDEM. Furthermore, disaggregating the results into percentages of tree cover highlights the improvements enabled by removing vegetation bias from COPDEM. For tree cover percentages ranging between 0% and 20%, COPDEM showed an RMSE of 6.40 m and a BIAS of 2.2 m, whereas the RMSE and BIAS of ANADEM were reduced to 5.83 m and 1.67 m, respectively. On the other hand, for areas with greater than 80% tree cover percentage, RMSE and BIAS were reduced from 17.20 and 15.38 m (COPDEM) to 6.94 and 0.25 m (ANADEM). More detailed results in terms of the metrics at different ranges of vegetation cover are presented in Table 2.

Previous studies have identified a significant error discrepancy in global DEMs within densely forested areas [47,48,49,50], primarily due to sensor limitations when accurately capturing terrain elevation. Li et al. [47] evaluated COPDEM in China, showing an RMSE of 11.9 m with a standard deviation of 8.52 m in regions predominantly covered by dense forests. Similarly, Hawker et al. [50] evaluated three models (MERIT, SRTM and TanDEM-X) with high accuracy LiDAR data, showing that all of these models were associated with high RMSE (10–14 meters).

Figure 3 shows the error distribution in the ANADEM and COPDEM models in relation to tree cover percentage. For areas with greater than 20% vegetation cover, the benefit of the bias removal performed to generate ANADEM becomes increasingly evident. At approximately 80% vegetation cover, COPDEM results are concentrated around a bias close to 20 m, whereas ANADEM results concentrate at values much closer to zero. This is consistent with the findings of Hawker et al. [20] for COPDEM in areas with canopy cover exceeding 50%, where an RMSE of 26 m was observed. On the other hand, ANADEM was found to be capable of decreasing error differences related to vegetation bias, resulting in a more centered normal distribution, despite exhibiting an increase in negative error differences.

3.2. Spatial Analysis

To evaluate the spatial results obtained by removing vegetation bias from COPDEM, we selected two areas to illustrate the improvements observed using ANADEM. Figure 4 shows a comparison between the original elevation data from COPDEM, data obtained ANADEM, and the extracted vegetation bias removed from COPDEM. Overall, ANADEM significantly reduced the influence of vegetation on the elevation values, resulting in a smoother surface than COPDEM. Although certain locations in ANADEM exhibited residual vegetation features due to an underestimation of the vegetation bias, the residual biases were considerably lower than those observed in COPDEM.

Similarly, Figure 5 provides a comparison of results for COPDEM, ANADEM, and the percentage of vegetation removed in an area of the Amazon with extensive land cover changes. The effective removal of vegetation in ANADEM is once again evident, leading to improvements in the representation of the terrain surface, although some vegetation features are still present. Uncertainties in vegetation removal in ANADEM may also occur due to seasonal variations in vegetation, which are influenced by the differences in the time periods of acquisition between COPDEM and multispectral data from satellites used to represent vegetation.

3.3. Evaluation of Drainage Network, Slope and Aspect Subsection

DEM inaccuracies can significantly influence the information derived, such as slope, aspect, and drainage network [51,52]. To evaluate the advances obtained with ANADEM due to vegetation removal, Figure 6 shows an example of a drainage network obtained from COPDEM and ANADEM. The drainage network was obtained using TerraHidro software provide by the Brazilian National Institute for Space Research (INPE) [53]. Initially, we removed topographic depressions and applied the eight-direction (D8) model to generate a flow direction map and a flow accumulation map, and then creating a drainage network by considering a threshold of 1000 cells.

Visual inspection demonstrates that ANADEM-delineated streams are more realistic and smoother than those delineated using COPDEM; this reflected the lower vegetations bias in ANADEM. Common errors in defining drainage networks using DEMs may result from the misinterpretation of flow direction due to objects on the surface acting as ‘barriers’, leading to an inaccurate delineation of the drainage network. Thus, the most significant differences were observed in areas with dense vegetation, which the drainage network computation failed to represent, as shown for COPDEM in Figure 6. This coincides with previous studies demonstrating a loss of drainage delineation accuracy in areas with high vegetation cover [54,55].

The vertical accuracy of both COPDEM and ANADEM is illustrated in Figure 7 in terms of slope (Figure 7a) and aspect (Figure 7b). The vertical error in slope, when compared with ICESat-2 data, showed a proportional increase with increasing slope for both models, although the improvements for ANADEM were more significant. In flat areas (slope <3°), ANADEM reduced the error by 84.9% (3.7 m to 0.57 m), whereas on steep slopes (>21°), ANADEM reduced the error by 44.5% (from 9.5 m to 5.3 m). For aspect, the errors were practically independent of slope orientation, with the lowest and highest errors observed in the W and NNE directions (5.97 m/6.48 m, respectively), whereas for ANADEM the S and E directions showed the lowest and highest errors (1.30 m/1.63 m, respectively). These results suggest that slope has a greater influence than aspect on the vertical accuracy of DEMs, supporting the conclusions of previous studies [56,57]. Hence, vegetation removal led to a reduction in errors for ANADEM without clearly indicating changes in dependency on slope and aspect.

3.4. Comparative Assessment with Globally Available Models

The accuracy of ANADEM was compared against terrain elevation from ICESat-2 and four others DEMs with free global availability. Table 3 shows the performance statistics for each model based on 600,000 randomly selected points. Among the models examined, SRTM and MERIT have the lowest spatial resolutions (90 meters), which influences their statistical performance. Canopy structure mainly influences SRTM performance, showing an RMSE of 16.86 m and a bias of 10.33 m. MERIT DEM shows a clear improvement in removing vegetation bias from the SRTM model, reducing RMSE and BIAS to 8.94 m and 3.88 m, respectively. Additionally, COPDEM, with a spatial resolution of 30 meters, showed an improvement in accuracy relative to SRTM. Despite the influence of the presence of vegetation, COPDEM results showed an average RMSE of 12.40 m and an average BIAS of 9.56 m.

The FABDEM and ANADEM models presented the lowest errors, denoting a clear improvement in terrain elevation representation compared to other models. The vegetation bias was lowest for ANADEM (1.50 m), whereas FABDEM showed an overestimation of 1.76 m. Our results determined an RMSE of 6.81 m for FABDEM, whereas ANADEM showed a slightly higher RMSE of 6.99 m. The median error compared to ICESat-2 data among the models ranged from 0.75 to 9.28 m. ANADEM exhibited the lowest median error. The average standard deviation error ranged from 6.47 m (FABDEM) to 12.75 m (SRTM).

Figure 8 characterizes the performance metrics into two groups of tree cover percentages. For points classified with less than 50% tree cover (Figure 8a), ANADEM exhibited the lowest error (1.30 m). The errors associated with FABDEM, MERIT, COPDEM, and SRTM were 1.62, 2.16, 2.79, and 4.27 m, respectively. As expected, higher errors were obtained for areas with greater than 50% tree cover. SRTM featured the highest errors (15.78 m), followed by COPDEM (14.85 m). MERIT decreased vegetation bias considerably (5.60 m) despite having a higher value than ANADEM and FABDEM. The most accurate results were obtained for ANADEM and FABDEM, with an average BIAS of 1.49 m and 2.16 m, respectively.

Validation results for each model were also categorized by land cover classes (Table 4) using the MapBiomas classification [58], which is available at a spatial resolution of 30 meters for Brazil. Results were obtained for most representative classes of the study area, totaling six classes: forests, urban areas, savannas, croplands, grasslands, and pastures. The highest bias was observed for forested areas, ranging from 0.4 m (ANADEM) to 14.31 m (COPDEM). ANADEM showed the lowest bias in pastures (1.42 m), grasslands (0.45 m), croplands (0.36 m), and forests, whereas FABDEM showed better results for savannas (2.17 m) and urban areas (0.17 m). The percentage RMSE for ANADEM varied from 2.26 m (croplands) to 7.09 m (forests), whereas those for COPDEM varied between 2.39 (croplands) to 16.44 m (forests). In urban areas, despite the lack of application of a specific bias correction caused by buildings in ANADEM, its performance metrics were slightly better than those of COPDEM. This may be attributed to the post-processing adopted in ANADEM to reduce noise generated by the machine learning model. Additionally, MERIT exhibited a lower performance relative to COPDEM in certain classifications, for example, croplands, grasslands, pastures, and urban areas.

4. Conclusions

In this study, we presented ANADEM, a DTM for South America based on a solution to remove the influence of vegetation on elevation values from COPDEM. ANADEM uses remote sensing data, cloud computation, and machine learning techniques. According to our findings, ANADEM reduced 85.6% of the vegetation effect relative to COPDEM, while also exhibiting accurate results in comparison with other models. Because vegetation influences stream delineation, additional benefits over COPDEM were observed when extracting the drainage network from ANADEM.

Several concluding points should be emphasizing. Firstly, the accuracy of ANADEM results is directly impacted by the terrain elevation values used for model training, as well as by the multispectral data representing vegetation. Although we utilized publicly available high-quality data in our study, errors may arise in areas with high cloud cover or inconsistent terrain elevation values from ICESat-2; these may introduce uncertainties in vegetation removal. Secondly, the parameters and temporal data period used to generate ANADEM were defined to facilitate vegetation removal on a continental scale, which could result in inconsistencies due to temporal differences in data acquisition. In such cases, users are encouraged to reproduce the proposed methodology by incorporating adaptations for smaller scales and by selecting multispectral data that align with the vegetation bias observed in COPDEM. Finally, we did not apply bias corrections for urban areas due to a perceived lack of publicly available datasets for model training that accurately represent buildings. Future studies should aim to better understand the biases introduced by the presence of buildings and explore methods of bias identification without losing terrain characteristics.

ANADEM data are freely available to users and we anticipate that the proposed methods may enable new advancements in hydrological and surface modelling in South America. Additionally, ANADEM may facilitate the development of new products based on the scripts utilized in this study for removing vegetation bias from DEM models.