AGB Estimation in a Tropical Mountain Forest (TMF) by Means of RGB and Multispectral Images Using an Unmanned Aerial Vehicle (UAV)

The present investigation evaluates the accuracy of estimating above-ground biomass (AGB) by means of two different sensors installed onboard an unmanned aerial vehicle (UAV) platform (DJI Inspire I) because the high costs of very high-resolution imagery provided by satellites or light detection and ranging (LiDAR) sensors often impede AGB estimation and the determination of other vegetation parameters. The sensors utilized included an RGB camera (ZENMUSE X3) and a multispectral camera (Parrot Sequoia), whose images were used for AGB estimation in a natural tropical mountain forest (TMF) in Southern Ecuador. The total area covered by the sensors included 80 ha at lower elevations characterized by a fast-changing topography and different vegetation covers. From the total area, a core study site of 24 ha was selected for AGB calculation, applying two different methods. The first method used the RGB images and applied the structure for motion (SfM) process to generate point clouds for a subsequent individual tree classification. Per the classification at tree level, tree height (H) and diameter at breast height (DBH) could be determined, which are necessary input parameters to calculate AGB (Mg ha−1) by means of a specific allometric equation for wet forests. The second method used the multispectral images to calculate the normalized difference vegetation index (NDVI), which is the basis for AGB estimation applying an equation for tropical evergreen forests. The obtained results were validated against a previous AGB estimation for the same area using LiDAR data. The study found two major results: (i) The NDVI-based AGB estimates obtained by multispectral drone imagery were less accurate due to the saturation effect in dense tropical forests, (ii) the photogrammetric approach using RGB images provided reliable AGB estimates comparable to expensive LiDAR surveys (R2: 0.85). However, the latter is only possible if an auxiliary digital terrain model (DTM) in very high resolution is available because in dense natural forests the terrain surface (DTM) is hardly detectable by passive sensors due to the canopy layer, which impedes ground detection.


Introduction
Accurate information about forest cover, land use and above-ground biomass (AGB) are critical parameters for many environmental studies as well as for conservation initiatives concerning the mitigation of global warming, such as REDD+ (Reducing emissions from deforestation and forest degradation and the role of conservation, sustainable management of forests and enhancement of forest carbon stocks in developing countries) [1,2]. By means of such data and their temporal development, The climate within the study area is per-humid with marked altitudinal gradients in air temperature, air humidity [57,58], cloudiness [59], rainfall [60], and wind conditions [56]. Mean annual air temperature ranges from 19.4 °C at the valley bottom to 9.4 °C at the mountain tops, whereas average relative humidity varies between 70% at the valley bottom and at open sites to nearly 100% at the mountain ridges and inside the TMF. Precipitation shows a clear annual cycle with a main rainy season in austral winter (between May and September) and a relative dry season in austral summer (between November and February). Wind direction is predominately from the east due to the tropical easterlies, reaching average monthly velocities up to 15.5 m s -1 at the ridges in austral winter [61]. The climate within the study area is per-humid with marked altitudinal gradients in air temperature, air humidity [57,58], cloudiness [59], rainfall [60], and wind conditions [56]. Mean annual air temperature ranges from 19.4 • C at the valley bottom to 9.4 • C at the mountain tops, whereas average relative humidity varies between 70% at the valley bottom and at open sites to nearly 100% at the mountain ridges and inside the TMF. Precipitation shows a clear annual cycle with a main rainy season in austral winter (between May and September) and a relative dry season in austral summer (between November and February). Wind direction is predominately from the east due to the tropical easterlies, reaching average monthly velocities up to 15.5 m s −1 at the ridges in austral winter [61].

Equipment
To estimate the AGB in the natural TMF, the UAV DJI Inspire 1 [48] was used, which was equipped with an RGB camera (ZENMUSE X3) and a multispectral camera (Parrot Sequoia) [62].

DJI Inspire 1
The UAV DJI Inspire 1 is a multi-rotor drone with 4 electric motors (quadcopter, Figure 2) that includes a remote controller that operates up to 2 km (radio) under unobstructed environmental and undisturbed meteorological conditions [48]. The maximum speed of the drone is 22 m s −1 , and its maximum operation altitude is 4500 m above sea level. Factory built, the DJI Inspire 1 is equipped with the RGB ZENMUSE X3 camera with a resolution of 12 megapixels, which also includes an integrated gimbal system that provides stability and avoids distortion and blurring. The flight time varies with respect to payload, altitude above sea level, and weather conditions (mainly wind speed) but generally lies between 12-18 min. By considering the additional weight of a multispectral sensor and its external battery (approximately 250 g), an average flight time of about 10 min could be achieved for the mountainous environment in the study area.
To estimate the AGB in the natural TMF, the UAV DJI Inspire 1 [48] was used, which was equipped with an RGB camera (ZENMUSE X3) and a multispectral camera (Parrot Sequoia) [62].

DJI Inspire 1
The UAV DJI Inspire 1 is a multi-rotor drone with 4 electric motors (quadcopter, Figure 2) that includes a remote controller that operates up to 2 km (radio) under unobstructed environmental and undisturbed meteorological conditions [48]. The maximum speed of the drone is 22 m s -1 , and its maximum operation altitude is 4500 m above sea level. Factory built, the DJI Inspire 1 is equipped with the RGB ZENMUSE X3 camera with a resolution of 12 megapixels, which also includes an integrated gimbal system that provides stability and avoids distortion and blurring. The flight time varies with respect to payload, altitude above sea level, and weather conditions (mainly wind speed) but generally lies between 12-18 min. By considering the additional weight of a multispectral sensor and its external battery (approximately 250 g), an average flight time of about 10 min could be achieved for the mountainous environment in the study area.

Parrot Sequoia
The Parrot Sequoia camera is a multispectral sensor [62] that allows for the capture of imagery at 4 specific bands within the visible and infrared electromagnetic spectrum: Green (550 nm), red (660 nm), red edge (735 nm), and near infrared (790 nm). Furthermore, the camera includes a sunshine module that automatically calibrates the received images to compensate for the variability in sunlight conditions during flight as well as during different campaigns [63]. The specific bands of the Parrot Sequoia camera are similar to multispectral sensors onboard Landsat and NOAA-AVHRR satellites, allowing for comparable investigations to be conducted [62].

Parrot Sequoia
The Parrot Sequoia camera is a multispectral sensor [62] that allows for the capture of imagery at 4 specific bands within the visible and infrared electromagnetic spectrum: Green (550 nm), red (660 nm), red edge (735 nm), and near infrared (790 nm). Furthermore, the camera includes a sunshine module that automatically calibrates the received images to compensate for the variability in sunlight conditions during flight as well as during different campaigns [63]. The specific bands of the Parrot Sequoia camera are similar to multispectral sensors onboard Landsat and NOAA-AVHRR satellites, allowing for comparable investigations to be conducted [62].

Flight Planning and Data Acquisition
Depending on the size of study area, specific characteristics of the employed UAV must be considered (especially payload and resulting flight time) to determine the area covered by one single flight as well as the number of flights necessary to cover the whole area of interest. For the present research, a core site of 24 ha was selected, located at lower elevations in the San Francisco catchment, including parts of the valley bottom (river), side valleys, and ridges covered by natural TMF in which small gaps of natural succession (ancient landslides, grassland, and scrubs) are also present ( Figure 3). To ensure complete coverage of the core site, the limits for data acquisition (vertical photographs) were extended at all borders (total area of coverage: 80 ha).
The flight characteristics were configured by means of the free software Precision Flight [64], an application that runs on an Android device that afterwards automatically executes the planned survey. The flight path was designed in an east-west direction with a nominal speed of 9 m s −1 to cover the whole area (80 ha) by two flights of approximately 10 min. The data were taken at a flight height of 300 m above the starting point with 90% side and forward overlap of the images, as recommended by the Pix4D documentation for locations with fast-changing topography and forests with dense canopies, to obtain enough points for individual image matching [65].
To improve the accuracy of the obtained images (coordinates X, Y and Z), ground control points (GCPs, Figure 3) were set within the area of coverage before executing the flights [28,66]. A minimum of 3 GCPs are necessary to ensure reliable image accuracy, but five to 10 GCPs are more suitable, particularly for larger areas [67]. It has been proven that a higher number of GCPs do not improve the final product significantly, nor do they improve the image accuracy [68], so for this reason six GCPs were set for the present study. The GCPs were located at flat open sites inside the total area (80 ha) to guarantee easy detection within the images. Therefore, dense canopies, steep slopes, and positions too close to the borders were avoided because the local topography as well as the camera position modifies the nadir angles of the different objects [29]. The exact geographical position of the GCPs ( Figure 3) were determined using a GPS Trimble R6 system. Furthermore, the multispectral camera needs a radiometric calibration target to guarantee reliable values for each individual spectral band. Therefore, an AirInov target, specifically designed for the Parrot Sequoia camera, was used to calibrate the bands before flight execution. More information about the target and its application can be found in [63]. Depending on the size of study area, specific characteristics of the employed UAV must be considered (especially payload and resulting flight time) to determine the area covered by one single flight as well as the number of flights necessary to cover the whole area of interest. For the present research, a core site of 24 ha was selected, located at lower elevations in the San Francisco catchment, including parts of the valley bottom (river), side valleys, and ridges covered by natural TMF in which small gaps of natural succession (ancient landslides, grassland, and scrubs) are also present ( Figure  3). To ensure complete coverage of the core site, the limits for data acquisition (vertical photographs) were extended at all borders (total area of coverage: 80 ha).
The flight characteristics were configured by means of the free software Precision Flight [64], an application that runs on an Android device that afterwards automatically executes the planned survey. The flight path was designed in an east-west direction with a nominal speed of 9 m s -1 to cover the whole area (80 ha) by two flights of approximately 10 min. The data were taken at a flight height of 300 m above the starting point with 90% side and forward overlap of the images, as recommended by the Pix4D documentation for locations with fast-changing topography and forests with dense canopies, to obtain enough points for individual image matching [65].
To improve the accuracy of the obtained images (coordinates X, Y and Z), ground control points (GCPs, Figure 3) were set within the area of coverage before executing the flights [28,66]. A minimum of 3 GCPs are necessary to ensure reliable image accuracy, but five to 10 GCPs are more suitable, particularly for larger areas [67]. It has been proven that a higher number of GCPs do not improve the final product significantly, nor do they improve the image accuracy [68], so for this reason six GCPs were set for the present study. The GCPs were located at flat open sites inside the total area (80 ha) to guarantee easy detection within the images. Therefore, dense canopies, steep slopes, and positions too close to the borders were avoided because the local topography as well as the camera position modifies the nadir angles of the different objects [29]. The exact geographical position of the GCPs ( Figure 3) were determined using a GPS Trimble R6 system. Furthermore, the multispectral camera needs a radiometric calibration target to guarantee reliable values for each individual spectral band. Therefore, an AirInov target, specifically designed for the Parrot Sequoia camera, was used to calibrate the bands before flight execution. More information about the target and its application can be found in [63]. The UAV flights were executed under sunny weather conditions during 26 April 2018, before the main rainy season started. During the survey, 132 RGB images and 230 multispectral images were obtained, which resulted in a total raw data size of 8 GB.

Data Processing
The data processing chain for AGB estimation by means of the RGB and multispectral images is shown in Figure 4. Both image types were analyzed independently. The UAV flights were executed under sunny weather conditions during 26 April 2018, before the main rainy season started. During the survey, 132 RGB images and 230 multispectral images were obtained, which resulted in a total raw data size of 8 GB.

Data Processing
The data processing chain for AGB estimation by means of the RGB and multispectral images is shown in Figure 4. Both image types were analyzed independently.

AGB Estimation by Means of RGB Data-Use of Photogrammetry
Allometric equations are widely used for AGB estimation [69][70][71] because they include forest structure parameters, such as tree height (H), tree diameter at breast height (DBH), and wood density (WD). The required information at tree level can be obtained by remote sensed data (H and DBH) and non-destructive field measurements (H, DBH, and WD), which do not harm the ecosystem [8,72]. For individual H detection, 3D data (point clouds) are required, which can be derived from expensive LiDAR surveys [73] but also from RGB data when applying computational photogrammetric techniques [1,74].
To obtain the 3D data, first, the raw data of the ZENMUSE X3 camera were checked and all images at the beginning and at the end of each individual flight eliminated because these images do not comply with the established flight characteristics (flight height, image overlap, and pixel resolution). The revised dataset included 132 images, which were used to generate an orthophoto ( Figure 3) and a point cloud for the whole study area ( Figure 5). The remaining images were processed using the Pix4D software [65], which uses the structure from motion (SfM) process [75]. This process applies the principles of traditional stereoscopic photogrammetry in which the

AGB Estimation by Means of RGB Data-Use of Photogrammetry
Allometric equations are widely used for AGB estimation [69][70][71] because they include forest structure parameters, such as tree height (H), tree diameter at breast height (DBH), and wood density (WD). The required information at tree level can be obtained by remote sensed data (H and DBH) and non-destructive field measurements (H, DBH, and WD), which do not harm the ecosystem [8,72]. For individual H detection, 3D data (point clouds) are required, which can be derived from expensive LiDAR surveys [73] but also from RGB data when applying computational photogrammetric techniques [1,74].
To obtain the 3D data, first, the raw data of the ZENMUSE X3 camera were checked and all images at the beginning and at the end of each individual flight eliminated because these images do not comply with the established flight characteristics (flight height, image overlap, and pixel resolution). The revised dataset included 132 images, which were used to generate an orthophoto ( Figure 3) and a point cloud for the whole study area ( Figure 5). The remaining images were processed using the Pix4D software [65], which uses the structure from motion (SfM) process [75]. This process applies the principles of traditional stereoscopic photogrammetry in which the horizontal and vertical position of specific geometrical features are determined using several images of the same object taken under different viewing angles [74]. This can be done by means of a series of RGB images with a side and forward overlap of between 60% and 90%. For the present study, an overlap of 90% was fixed to avoid mismatches during the overlap process of the images. Furthermore, the fast-changing topography and forest structure in the study area requires a higher overlap to obtain enough points for image matching [65].  As Aasen et al. [28] stated, the SfM processing provides a high geometric fidelity if the obtained/captured images contain the exact geographical position and orientation of the spectral sensor. These conditions are given for the present study because sensor orientation is implemented in the ZENMUSE X3 camera (orientation), and the DJI Inspire 1 drone includes an on-board inertial measurement unit (IMU) sensor receiving GNSS signals (geographical position) [48]. As mentioned before, to improve the accuracy of the resulting images, six GCPs were set within the study area, which guarantees reliable image accuracy [76]. The location of the GCPs was introduced in the photogrammetric process, which identifies the specific GCP in each individual image and merges the image accordingly. The six GCPs were identified in 120 images, which is why a high geometric fidelity of the resulting information can be assumed [28,67].
In general, the workflow of the SfM approach consists of three major steps: (i) Feature extraction, (ii) feature matching, and (iii) reconstruction. Most of the processing time is needed for feature matching [76]. The software overlaps all available images and generates 3D-point clouds by means of the different angle views of the individual images [65,76]. Therefore, sufficient GCPs are necessary, especially in primary forests, to obtain reliable image matching [66]. The obtained point cloud was reduced to the core study site, located in the center of the investigation area. The Pix4D software also generates an orthophoto whose resolution depends on flight altitude (here: 300 m above ground). The programmed flight characteristics resulted in a pixel size (GSD) of about 25 cm for the orthophoto (see Figure 3) and the 3D point cloud ( Figure 5).
To obtain H and DBH estimates for individual trees, the 3D point cloud was processed according to the method presented by González-Jaramillo et al. [8]. From the 3D point cloud, a digital terrain model (DTM) and a digital surface model (DSM) were created, in which the DTM model can be interchanged with any other DTM in very high resolution available for the study area [77]. The utilized software (Fusion 3.7) [78] subtracts the DTM from the DSM to obtain a canopy height model (CHM), which displays the difference between the ground (terrain, bare-earth points) and the highest elevation returns for each GSD [79]. In this case, the CHM depicts the vegetation height while the Canopy Maxima method was also applied within a variable-size evaluation window (here: 3 × 3) [80] to estimate individual H of dominant trees, which represent 70% to 90% of the total AGB of a forest As Aasen et al. [28] stated, the SfM processing provides a high geometric fidelity if the obtained/captured images contain the exact geographical position and orientation of the spectral sensor. These conditions are given for the present study because sensor orientation is implemented in the ZENMUSE X3 camera (orientation), and the DJI Inspire 1 drone includes an on-board inertial measurement unit (IMU) sensor receiving GNSS signals (geographical position) [48]. As mentioned before, to improve the accuracy of the resulting images, six GCPs were set within the study area, which guarantees reliable image accuracy [76]. The location of the GCPs was introduced in the photogrammetric process, which identifies the specific GCP in each individual image and merges the image accordingly. The six GCPs were identified in 120 images, which is why a high geometric fidelity of the resulting information can be assumed [28,67].
In general, the workflow of the SfM approach consists of three major steps: (i) Feature extraction, (ii) feature matching, and (iii) reconstruction. Most of the processing time is needed for feature matching [76]. The software overlaps all available images and generates 3D-point clouds by means of the different angle views of the individual images [65,76]. Therefore, sufficient GCPs are necessary, especially in primary forests, to obtain reliable image matching [66]. The obtained point cloud was reduced to the core study site, located in the center of the investigation area. The Pix4D software also generates an orthophoto whose resolution depends on flight altitude (here: 300 m above ground). The programmed flight characteristics resulted in a pixel size (GSD) of about 25 cm for the orthophoto (see Figure 3) and the 3D point cloud ( Figure 5).
To obtain H and DBH estimates for individual trees, the 3D point cloud was processed according to the method presented by González-Jaramillo et al. [8]. From the 3D point cloud, a digital terrain model (DTM) and a digital surface model (DSM) were created, in which the DTM model can be interchanged with any other DTM in very high resolution available for the study area [77]. The utilized software (Fusion 3.7) [78] subtracts the DTM from the DSM to obtain a canopy height model (CHM), which displays the difference between the ground (terrain, bare-earth points) and the highest elevation returns for each GSD [79]. In this case, the CHM depicts the vegetation height while the Canopy Maxima method was also applied within a variable-size evaluation window (here: 3 × 3) [80] to estimate individual H of dominant trees, which represent 70% to 90% of the total AGB of a forest stand [18,19]. For this study, only dominant trees higher than 5 m with a DBH greater than 10 cm were considered to calculate the AGB, as Gianico et al. [79] recommended.
The H of the individual trees detected was used to calculate the DBH based on the height-diameter relationship equation established by González-Jaramillo et al. [8] for the same study area (Equation (1)), where DBH is the diameter at breast height in cm, and H is the estimated height of the tree (m), obtained from the point cloud.
To estimate the AGB (Mg ha −1 ), mean WD from a previous investigation was used, which was obtained by means of field measurements taken in the study area (0.59 g cm −3 ) [81]. The applied AGB equation at tree level corresponds to Chave et al. [49], who suggested a specific equation for tropical wet forest, which can be written as follows (Equation (2)), where AGB tree corresponds to the AGB of a specific tree (Mg), WD represents the wood density average (gr cm −3 ), H is the obtained height of each detected tree (m) and D is the estimated DBH (cm) obtained by means of Equation (1). The product applying Equation (2) is the AGB per individual tree. Finally, a grid layer with a 1 ha × 1 ha resolution was overlaid to determine AGB per hectare, where all individual tree AGB was added up for each grid cell [8].

AGB Estimation Using Multispectral Data
To estimate AGB by means of the multispectral data obtained from the Parrot Sequoia camera, the equation proposed by Das and Singh [50] was applied. They presented a specific equation for different tropical forest types. The selected equation corresponds to tropical evergreen forest, which is comparable to the forest type presented in the study area. The raw data with a spatial resolution of 25 cm (multispectral imagery) were processed using the Pix4D software [65] which allows radiometric and geometric corrections of each spectral band. As mentioned before, the radiometric calibration was done using an AirInov target [63] before executing the flights. The specific albedo values were provided by the manufacturer [62], which allowed for the correction of the reflectivity values, ranging between 0 and 100.
To distinguish between forest, bare soil, and shrubs in the corrected images, first, a false color composite was generated. Therefore, the spectral bands NIR, RED, and GREEN were used [82], which highlight areas with vegetation in red, whereas other ground covers are shown in a different color. The false color composite was the basis for a non-supervised classification, using the Iso Cluster Unsupervised Classification tool available in ArcGis 10.5.1. [83,84]. This tool evaluates the whole dataset and classifies the vegetation into different categories (here: Herbs/bare soil, scrubland, and forest).
Then, the normalized difference vegetation index (NDVI) [85] was calculated, which expresses the vigorousness of the vegetation. The NDVI is directly related to the photosynthetic capacity and Remote Sens. 2019, 11, 1413 9 of 22 therefore to the energy absorption of the vegetation. Its magnitude ranges between −1 and 1, in which negative values indicate water bodies, values near 0 bare soils, and positive values vegetation cover [86]. The NDVI is determined by the ratio of the NIR band (near infrared, 790 nm) and the RED band (red, 660 nm), expressed in Equation (3), where NDVI is the normalized difference vegetation index, NIR = near-infrared band, and RED = red spectral band. After the NDVI calculation, a regular grid mask of 1 hectare was overlaid, and all 25 cm pixels within a grid averaged to obtain a mean NDVI value per hectare. Finally, the equation for evergreen forest proposed by Das and Singh [50] was applied to calculate the AGB (Mg ha −1 ), which can be written as follows, AGB = 324.2 × NDVI mean + 14.18 (4) where AGB is the above-ground biomass (Mg ha −1 ), and NDVI mean is the normalized difference vegetation index, averaged for 1 ha.

Validation of the RGB and Multispectral AGB Estimations
The two different techniques for AGB estimation (RGB and multispectral data) using a UAV were validated by means of an independent AGB dataset (Mg ha −1 ) derived from LiDAR data. The LiDAR survey was executed in the same study area in 2012, using a Leica Geosystem ALS-50-II CM laser scanner installed onboard a Eurocopter AS350B2 Ecuriel Helicopter. The resulting point cloud density was at least 10 pulses per 1 m 2 [31], from which a DTM, DSM, and CHM model in a 25 cm x 25 cm resolution was generated, and the AGB (Mg ha −1 ) calculated applying an allometric equation [49]. For more details about the LiDAR-AGB calculation, please refer to González-Jaramillo et al. [8].
The RGB model executed in this study is similar to the LiDAR approach [8], although the RGB camera is a passive sensor, whereas the LiDAR laser scanner an active sensor. Therefore, the generated DTM, DSM, CHM and individual tree parameters obtained from the RGB data were compared to the reference data from LiDAR. The multispectral model is different because the AGB (Mg ha −1 ) is directly calculated by means of the equation for tropical evergreen forests proposed by Das and Singh [50], which is why only the resultant AGB values were compared to the reference LiDAR-AGB.
The accuracy of the obtained results, RGB or multispectral vs. LiDAR, were determined by means of the coefficient of determination (R 2 ) and the root mean square error (RMSE) [46,87]. The equation for R 2 (Equation (5)) and RMSE (Equation (6)) are written as follows, where x i and y i are the estimated and measured values, x and y are the average estimated and measured values, and n is the total number of existing values with respect to the compared parameters.

AGB Results by Means of RGB Data
By executing the explained photogrammetric processes for the RGB data, an orthophoto ( Figure 3) and a point cloud image were obtained ( Figure 5). The orthophoto was used simply for supervision, while from the point cloud, the horizontal and vertical structure of the TMF at tree level was determined.
The point cloud image had an average density of 30 points/m 2 for the total zone (~80 ha), but zones with no information were also present, especially near the borders of the flight domain, for which reason a core area (24 ha) was established to avoid gaps in the information. Furthermore, due to the dense canopies and the fast-changing topography in the core area, which impede the determination of bare earth points by means of images of the RGB information, the high resolution DTM (25 cm × 25 cm) generated from the LiDAR survey (2012) was used to generate the CHM. As Karpina et al. [77] stated, this interchange is practicable because passive remote sensors cannot detect the terrain surface under dense canopy layers. Therefore, to generate the CHM, the ancillary LiDAR-DTM was subtracted from the RGB-DSM. Then, the Canopy Maxima tool of FUSION 3.7 was executed to detect the local maxima in the CHM, which represents the H of individual trees and their location. A total of 7075 dominant trees were detected (Figure 6a), and, by means of Equation (1), their individual DBH calculated. To estimate AGB (Mg ha −1 ), a grid mask of 1 ha × 1 ha was overlaid and all AGB values of the individual trees within a grid cell were added up. The resulting AGB map of the core area is shown in Figure 6b.
their individual DBH calculated. To estimate AGB (Mg ha -1 ), a grid mask of 1 ha × 1 ha was overlaid and all AGB values of the individual trees within a grid cell were added up. The resulting AGB map of the core area is shown in Figure 6b.
H varied between 8.50 m to 35.32 m in the core area in which the individual tree distribution depended on the topographical position. Bigger trees were found in depressions (ravine forest) because these locations are topographically more protected compared to the ridges (ridge forest), and trees are generally taller there [22,88]. This is also caused by the prevailing soil types present at the different locations [53]. In general, the down-slope fluxes accumulate material and nutrients in the depression [22], which is why soil depth and nutrient concentrations are generally higher there, supporting the resulting tree growth [22,88]. The ridge forest parts are frequently affected by soil erosion processes due to the steeper slopes and harsher climate conditions [56], which lead to shallow soils and unstable conditions (landslides) at these areas [55], subsequently reducing tree growth.
The DBH ranged between 10.01 cm and 127.39 cm, while distribution also depended on the topographical positions [54]. In general, the estimated AGB (Mg ha -1 ) was related to the topographical position and the specific forest type (ravine or ridge) and ranged between 18.77 Mg ha -1 (open ridge forest) and 317.77 Mg ha -1 (dense ravine forest), with a mean value of 148.83 Mg ha -1 for the core study area. In Figure 6, the aforementioned distribution is clearly visible, where higher AGB values (green and yellow colors) are located inside the side valley and near the valley bottom, whereas lower AGB values (orange and red colors) are displayed at the ridges or steep slopes, where small vegetation (shrubs or small trees) or regeneration areas prevail.  H varied between 8.50 m to 35.32 m in the core area in which the individual tree distribution depended on the topographical position. Bigger trees were found in depressions (ravine forest) because these locations are topographically more protected compared to the ridges (ridge forest), and trees are generally taller there [22,88]. This is also caused by the prevailing soil types present at the different locations [53]. In general, the down-slope fluxes accumulate material and nutrients in the depression [22], which is why soil depth and nutrient concentrations are generally higher there, supporting the resulting tree growth [22,88]. The ridge forest parts are frequently affected by soil erosion processes due to the steeper slopes and harsher climate conditions [56], which lead to shallow soils and unstable conditions (landslides) at these areas [55], subsequently reducing tree growth.
The DBH ranged between 10.01 cm and 127.39 cm, while distribution also depended on the topographical positions [54]. In general, the estimated AGB (Mg ha −1 ) was related to the topographical position and the specific forest type (ravine or ridge) and ranged between 18.77 Mg ha −1 (open ridge forest) and 317.77 Mg ha −1 (dense ravine forest), with a mean value of 148.83 Mg ha −1 for the core study area. In Figure 6, the aforementioned distribution is clearly visible, where higher AGB values (green and yellow colors) are located inside the side valley and near the valley bottom, whereas lower AGB values (orange and red colors) are displayed at the ridges or steep slopes, where small vegetation (shrubs or small trees) or regeneration areas prevail.

AGB Results by Means of Multispectral Data
By utilizing the multispectral data, a false color composite was generated as well as the individual NDVI values for each pixel (25 cm) calculated. The false color composite provided the basis for an unsupervised classification to divide the vegetation cover in herbs/bare soil, scrubland, and forest ( Figure 7a). The vegetation classification was used for visual inspection to verify the calculated NDVI values corresponding to these vegetation classes. The calculated NDVI map is illustrated in Figure 7b, displaying generally high values of between 0.38 and 0.86. Only water bodies (river course) showed negative values, while bare soil's NDVIs up to 0.20. The remaining land cover units had high positive values, which indicates vigorous vegetation [40,89]. Nonetheless, the forest distribution is detectable in the NDVI map, where higher values (dark green) were displayed inside the side valleys and slightly lower values (light green) at the ridges (Figure 7b). The individual NDVI values were averaged for each hectare and integrated in Equation (4)

AGB Results by Means of Multispectral Data
By utilizing the multispectral data, a false color composite was generated as well as the individual NDVI values for each pixel (25 cm) calculated. The false color composite provided the basis for an unsupervised classification to divide the vegetation cover in herbs/bare soil, scrubland, and forest (Figure 7a). The vegetation classification was used for visual inspection to verify the calculated NDVI values corresponding to these vegetation classes. The calculated NDVI map is illustrated in Figure 7b, displaying generally high values of between 0.38 and 0.86. Only water bodies (river course) showed negative values, while bare soil's NDVIs up to 0.20. The remaining land cover units had high positive values, which indicates vigorous vegetation [40,89]. Nonetheless, the forest distribution is detectable in the NDVI map, where higher values (dark green) were displayed inside the side valleys and slightly lower values (light green) at the ridges (Figure 7b). The individual NDVI values were averaged for each hectare and integrated in Equation (4) to obtain AGB (Mg ha -1 ). The resulting AGB map is illustrated in Figure 8, indicating values between 191.46 MG ha -1 and 252.11 Mg ha -1 . The mean value corresponds to 237.21 Mg ha -1 . However, due to the generally high NDVI in the study area (natural TMF), the differences between ravine and ridge forest were indistinct, which resulted in high AGB values for the whole area. The resulting AGB map is illustrated in Figure 8, indicating values between 191.46 MG ha −1 and 252.11 Mg ha −1 . The mean value corresponds to 237.21 Mg ha −1 . However, due to the generally high NDVI in the study area (natural TMF), the differences between ravine and ridge forest were indistinct, which resulted in high AGB values for the whole area.

Validation
The validation was realized by means of the statistic software R-Studio (version 1.1.453) [90]. First, the results obtained from the RGB camera were analyzed, in which the individual generated models (DTM, DSM, and CHM) were compared to the LiDAR models [8]. The DTM comparison was made for two transects (Figure 3, yellow lines), one over dense vegetation cover and another over sparse or non-vegetation sites (Figure 9, red and black lines). It is clearly visible in Figure 9a that the RGB camera can hardly detect the terrain surface under dense vegetation (differences up to 20 m) because this passive sensor is not able to penetrate dense canopy layers [91,92]. In contrast, for sparse or non-vegetated areas (Figure 9b), only small differences between the RGB-DTM and the LiDAR-DTM were obtained (differences up to 2 m). Therefore, the RGB-DTM generation is largely inaccurate for natural tropical forests, especially for areas with dense vegetation, where high canopies obstruct the detection of the terrain surface (RGB-DTM RMSE: ~9 m, Table 1).
In contrast, the DSM comparison for both transects showed good accordance between the RGB and the LiDAR model because the top of the surface (here: Crown level) can be analyzed by means of RGB data (Figure 9, green and grey lines). This is confirmed by the statistics (Table 1), which determined a coefficient of determination (R 2 ) of 0.99 with an RMSE of 3.05 m at very high significance level <0.001 (p-value). However, differences between dense vegetation (Figure 9a, R 2 of 0.98 and a RMSE of 1.38 m) and less-vegetated sites (Figure 9b, R 2 of 0.99 with an RMSE of 0.49 m) were notable, which is due to the irregular canopy layer in dense natural TMF [22] and the time span between the LiDAR and UAV survey (vegetation growth).

Validation
The validation was realized by means of the statistic software R-Studio (version 1.1.453) [90]. First, the results obtained from the RGB camera were analyzed, in which the individual generated models (DTM, DSM, and CHM) were compared to the LiDAR models [8]. The DTM comparison was made for two transects (Figure 3, yellow lines), one over dense vegetation cover and another over sparse or non-vegetation sites (Figure 9, red and black lines). It is clearly visible in Figure 9a that the RGB camera can hardly detect the terrain surface under dense vegetation (differences up to 20 m) because this passive sensor is not able to penetrate dense canopy layers [91,92]. In contrast, for sparse or non-vegetated areas (Figure 9b), only small differences between the RGB-DTM and the LiDAR-DTM were obtained (differences up to 2 m). Therefore, the RGB-DTM generation is largely inaccurate for natural tropical forests, especially for areas with dense vegetation, where high canopies obstruct the detection of the terrain surface (RGB-DTM RMSE:~9 m, Table 1).
In contrast, the DSM comparison for both transects showed good accordance between the RGB and the LiDAR model because the top of the surface (here: Crown level) can be analyzed by means of RGB data (Figure 9, green and grey lines). This is confirmed by the statistics (Table 1), which determined a coefficient of determination (R 2 ) of 0.99 with an RMSE of 3.05 m at very high significance level <0.001 (p-value). However, differences between dense vegetation (Figure 9a, R 2 of 0.98 and a RMSE of 1.38 m) and less-vegetated sites (Figure 9b, R 2 of 0.99 with an RMSE of 0.49 m) were notable, which is due to the irregular canopy layer in dense natural TMF [22] and the time span between the LiDAR and UAV survey (vegetation growth). Remote Sens. 2019, 11, x FOR PEER REVIEW 13 of 22 By means of the DTM and the DSM, a CHM could be generated, but due to the high errors in the RGB-DTM, the RGB-CHM error also showed high inaccuracy (RMSE = 8.65 m), and no correlation between RGB-CHM and LiDAR-CHM could be determined (R 2 = 0.18, Table 1). Therefore, as Karpina et al. [77] suggested, the RGB-DTM was replaced by an accurate DTM in very high resolution (LiDAR-DTM) to guarantee a reliable CHM generation and a subsequent individual tree classification. Using the LiDAR-DTM in combination with the RGB-DSM (RGB-CHM*, Table 1), errors were notably reduced (RMSE = 3.00 m) and a good correlation between the models was obtained (R 2 = 0.80). Therefore, the RGB-CHM* model was used for the individual tree classification to determine H and DBH for the final AGB estimation. Within the core study area, a total of 7075 dominant trees were detected on basis of the RGB-CHM*, applying the canopy maxima method [80]. In comparison to the LiDAR-CHM in which 7317 trees were detected for the same area, results were similar (96.69%, Table 2). From the individual tree detection, H was deviated and DBH of each tree calculated. As shown in Table 2, H and DBH values for the individual trees were in the same range comparing RGB and LiDAR data, which resulted in similar AGB values. However, maximum AGB values were slightly higher for the RGB data, which may be due to the RMSE of the RGB-DSM (Table 1) or the time span between the surveys (vegetation growth). By means of the DTM and the DSM, a CHM could be generated, but due to the high errors in the RGB-DTM, the RGB-CHM error also showed high inaccuracy (RMSE = 8.65 m), and no correlation between RGB-CHM and LiDAR-CHM could be determined (R 2 = 0.18, Table 1). Therefore, as Karpina et al. [77] suggested, the RGB-DTM was replaced by an accurate DTM in very high resolution (LiDAR-DTM) to guarantee a reliable CHM generation and a subsequent individual tree classification. Using the LiDAR-DTM in combination with the RGB-DSM (RGB-CHM*, Table 1), errors were notably reduced (RMSE = 3.00 m) and a good correlation between the models was obtained (R 2 = 0.80). Therefore, the RGB-CHM* model was used for the individual tree classification to determine H and DBH for the final AGB estimation. Within the core study area, a total of 7075 dominant trees were detected on basis of the RGB-CHM*, applying the canopy maxima method [80]. In comparison to the LiDAR-CHM in which 7317 trees were detected for the same area, results were similar (96.69%, Table 2). From the individual tree detection, H was deviated and DBH of each tree calculated. As shown in Table 2, H and DBH values for the individual trees were in the same range comparing RGB and LiDAR data, which resulted in similar AGB values. However, maximum AGB values were slightly higher for the RGB data, which may be due to the RMSE of the RGB-DSM (Table 1) or the time span between the surveys (vegetation growth). Furthermore, the RGB results (H, DBH, and AGB) were plotted against the LiDAR results ( Figure 10). Therefore, the deviated H ( Figure 10a) and calculated DBH values of all detected trees from the RGB (7075 individuals) were used and a linear regression analysis executed. The results showed a good correlation for H (R 2 = 0.83) and DBH (R 2 = 0.83). AGB amounts of the specific hectares (24 ha) were also compared, which resulted in a good correlation between RGB-AGB and LiDAR-AGB (R 2 = 0.85, at a significance level p-value < 0.001, Figure 10b).  Furthermore, the RGB results (H, DBH, and AGB) were plotted against the LiDAR results ( Figure 10). Therefore, the deviated H ( Figure 10a) and calculated DBH values of all detected trees from the RGB (7075 individuals) were used and a linear regression analysis executed. The results showed a good correlation for H (R² = 0.83) and DBH (R² = 0.83). AGB amounts of the specific hectares (24 ha) were also compared, which resulted in a good correlation between RGB-AGB and LiDAR-AGB (R² = 0.85, at a significance level p-value < 0.001, Figure 10b). Then, the AGB results obtained from the multispectral data were validated against the LiDAR estimates. In Table 3, the calculated minimum, mean, and maximum AGB values are presented, which indicate notable differences. The calculated minimum AGB of the multispectral data (191.46 Mg ha -1 ) was notably higher, whereas as maximum AGB (252.11 Mg ha -1 ) was lower than the LiDAR estimates. The absence of correlation between the multispectral and LiDAR results (p-value > 0.05) was confirmed by the very low R 2 (p-value < 0.01) and the very high RMSE (127.35 m). The large differences were due to the tendency towards saturation of multispectral sensors for healthy and vigorous vegetation [82], for which reason, generally, high NDVI values were obtained for the whole area (Figure 7b), which resulted in similar AGB values for all hectares under study and consequently in unclear forest structure detection.

Discussion
Regarding the UAV technology used for AGB estimation in this natural TMF in Southern Ecuador, two main issues could be identified. First, flight time and image resolution are related to topographic conditions as well as to the payload capacity of the employed UAV, which determine the area covered by one single flight as well as the final GSD resolution (pixel size). Furthermore, a high overlap of the individual images is necessary and an adequate number of GCPs must be established [67,68] to obtain the required image accuracy during the photogrammetric process. The Then, the AGB results obtained from the multispectral data were validated against the LiDAR estimates. In Table 3, the calculated minimum, mean, and maximum AGB values are presented, which indicate notable differences. The calculated minimum AGB of the multispectral data (191.46 Mg ha −1 ) was notably higher, whereas as maximum AGB (252.11 Mg ha −1 ) was lower than the LiDAR estimates. The absence of correlation between the multispectral and LiDAR results (p-value > 0.05) was confirmed by the very low R 2 (p-value < 0.01) and the very high RMSE (127.35 m). The large differences were due to the tendency towards saturation of multispectral sensors for healthy and vigorous vegetation [82], for which reason, generally, high NDVI values were obtained for the whole area (Figure 7b), which resulted in similar AGB values for all hectares under study and consequently in unclear forest structure detection.

Discussion
Regarding the UAV technology used for AGB estimation in this natural TMF in Southern Ecuador, two main issues could be identified. First, flight time and image resolution are related to topographic conditions as well as to the payload capacity of the employed UAV, which determine the area covered by one single flight as well as the final GSD resolution (pixel size). Furthermore, a high overlap of the individual images is necessary and an adequate number of GCPs must be established [67,68] to obtain the required image accuracy during the photogrammetric process. The second issue refers particularly to the RGB images and the 3D point cloud generation. As the study showed, terrain or ground data (DTM) under dense forest canopies cannot be detected accurately because passive sensors do not have the capacity to penetrate dense forest canopies [91,92]. However, the methodology can be applied in managed forest, where gaps are big enough and generally homogeneous, which permits the detection of the terrain [93]. In dense natural forest, an adequate auxiliary DTM in high resolution (here: LiDAR-DTM) is necessary to generate a realistic CHM for further analysis [77,94].
Considering the second issue, it could be stated that from the RGB-CHM* the structural parameters of the forest could be identified, which include the forest type (here: Ravine or ridge) as well as H of each individual tree detected. Biggest trees were found in the protected ravine TMF parts (valleys and side valleys) with H up to 35.32 m, similar to González-Jaramillo et al. [8], who used LiDAR data for individual tree classification (36.31 m, Table 1). The same maximum H was reported by Paulick et al. [54] and Leuschner et al. [95] within the same study area. Besides this, mean H (14.66 m) was equal to the values presented by of Homeier et al. [96], who installed field plots at lower elevation, which confirms the accuracy of the applied method.
By means of the individual H and the specific equation proposed by González-Jaramillo et al. [8] for the San Francisco catchment, DBH was calculated at tree level. Due to similar H values (Figure 10a), a good correlation for DBH between RGB and LiDAR data could be recorded (R 2 : 0.83, p-value < 0.001, Figure 10b). Mean DBH was only slightly higher (28.08 cm) than the validated values from the LiDAR data (27.57 cm, Table 1). However, other studies presented notably lower mean DBH for this TMF (Mean DBH: 9.8 cm) [97], but these investigations included mid-and understory trees as well as trees at higher elevations, where harsher climate conditions and the more exposed position reduce H and DBH of the trees [55,56].
The individual tree classification applied here only detects dominant trees because smaller mid-and understory trees are often not visible due to the upper canopy layer [78]. Therefore, RGB information, but also LiDAR approaches, slightly overestimate the real mean DBH in dense natural tropical forests [8]. Nonetheless, dominant trees are most important for AGB estimation because they represent 70% to 90% of the total AGB of a forest stand [18,19]. To calculate the AGB (Mg ha −1 ) by means of the individual tree classification, the specific allometric equation for wet forests from Chave et al. [49] was applied. The results clearly reflected the forest structure of this TMF because highest AGB was estimated inside the side valleys (ravine forest, up to 317.77 Mg ha −1 ), whereas lowest AGB was near the ridges (ridge forest, 18.77 Mg ha −1 , Figure 6). Mean AGB for the core area (148.66 Mg ha −1 , Table 1) was in the same range as the LiDAR approach [8], but also similar to other AGB estimations published in previous studies (150.0 Mg ha −1 ) based on field plot measurements at lower elevation [22,54,95,98]. The small AGB differences between RGB-AGB and LiDAR-AGB were due to the RMSE of the RGB-CHM* and the number of trees detected, considering also the time span between RGB and LiDAR survey (see Section 4.3). Also of note, the TMF in the study area is protected and consists of mature forest (60 to 80 years old trees), where trees generally grow slowly, especially the older individuals [88]. Therefore, bigger AGB increments over time could not be expected, but other land cover changes, like landslides, which occur naturally, provoke a reduction in AGB for the affected areas until natural regeneration reestablishes the vegetation cover.
The multispectral image approach to calculate AGB by means of the deviated NDVI and applying the equation from Das and Singh [50] resulted in marked overestimations, especially for less vegetated sites (Table 3). This was principally due to the tendency towards saturation of the spectral bands over dense forest covers [82], resulting in similar and high AGB values for all hectares, which is why forest structure could not be depicted. This type of sensor saturation was also reported by Gu et al. [99], who obtained very high NDVI values in areas with dense canopies in the USA. Therefore, the NDVI can hardly be used in natural TMF for AGB estimation because differences in the fast-changing topography and vegetation cover cannot be detected. The applied equation from Das and Singh [50], who established their NDVI-AGB relationship in an evergreen forest in the Western Ghat region of India, also might not be suitable for AGB estimation in this TMF in Southern Ecuador. Therefore, the equation should be adjusted to provide explicit and more accurate NDVI-AGB relationships for individual ecosystems, but the matter with spectral band saturation for dense forest stand will persist. Nonetheless, the multispectral approach might be useful for crops in combination with crop height information to calculate the AGB [46,100] due to the generally more uniform ground cover and topography.
In summary, the more effective and accurate approach to estimate AGB by means of UAVs in natural TMF is the RGB alternative, applying an analysis at tree level. During the process, forest structure as well as individual H and DBH of each tree detected can be determined [101]. These parameters, in combination with mean WD of the specific forest type, allowed a realistic AGB calculation. The RGB data permit the precise detection of the surface (DSM), but ancillary terrain data in high resolution (DTM) from other sources is necessary in natural tropical forest because the dense canopies impede ground detection (canopy closure often 100%) [39,58]. However, in natural mid and high latitude forests [102] as well as in managed forests [93], the RGB approach might be a cost-effective alternative to expensive LiDAR surveys [103] due to the different forest stand characteristics with open canopies or regular distance between individual trees, which allows for accurate ground detection [77].

Conclusions
UAVs have the potential to generate terrain and surface information in high temporal and spatial resolution, due to their portability and flexibility. Furthermore, UAVs can be equipped with different sensors, which allow a wide range of applications and contribute to the advances in the remote sensing field. In contrast to classical remote sensing data (satellite images), UAVs also avoid the problems of cloudiness, which is particularly advantageous in tropical high mountains. Nonetheless, their range is limited due to the battery capacity and the additional payload, which reduce their flight time. Therefore, UAVs can adequately monitor smaller areas, besides the ability to fill gaps in existing imagery (e.g., LiDAR data).
As the study showed, for a reliable AGB estimation (Mg ha −1 ) in natural TMF, a tree level classification is necessary. This can be reached by RGB images with high side and forward overlap (80% or higher), to obtain multiple viewing angles of the objects to generate 3D point clouds. However, in natural tropical forests, terrain information is difficult to capture with this passive sensor due to the dense canopy layer, which impedes ground detection. Therefore, auxiliary DTMs in high resolution must be considered to generate reliable CHMs for the subsequent classification at tree level and AGB calculation.
In contrast, multispectral images suffer saturation of the spectral bands over dense natural forest stands, which results in generally high NDVI values with small differences between the land cover units. Therefore, AGB estimation, based on NDVI values, overestimates the real amounts in natural TMF, especially at less vegetated sites.