Developing Allometric Equations for Teak Plantations Located in the Coastal Region of Ecuador from Terrestrial Laser Scanning Data

: Traditional studies aimed at developing allometric models to estimate dry above-ground biomass (AGB) and other tree-level variables, such as tree stem commercial volume (TSCV) or tree stem volume (TSV), usually involves cutting down the trees. Although this method has low uncertainty, it is quite costly and inefficient since it requires a very time-consuming field work. In order to assist in data collection and processing, remote sensing is allowing the application of nondestructive sampling methods such as that based on terrestrial laser scanning (TLS). In this work, TLS-derived point clouds were used to digitally reconstruct the tree stem of a set of teak trees ( Tectona grandis Linn. F.) from 58 circular reference plots of 18 m radius belonging to three different plantations located in the Coastal Region of Ecuador. After manually selecting the appropriate trees from the entire sample, semi-automatic data processing was performed to provide measurements of TSCV and TSV, together with estimates of AGB values at tree level. These observed values were used to develop allometric models, based on diameter at breast height (DBH), total tree height (h), or the metric DBH 2 × h, by applying a robust regression method to remove likely outliers. Results showed that the developed allometric models performed reasonably well, especially those based on the metric DBH 2 × h, providing low bias estimates and relative RMSE values of 21.60% and 16.41% for TSCV and TSV, respectively. Allometric models only based on tree height were derived from replacing DBH by h in the expression DBH 2 x h, according to adjusted expressions depending on DBH classes (ranges of DBH). This finding can facilitate the obtaining of variables such as AGB (carbon stock) and commercial volume of wood over teak plantations in the Coastal Region of Ecuador from only knowing the tree height, constituting a promising method to address large-scale teak plantations monitoring from the canopy height models derived from digital aerial stereophotogrammetry.


Introduction
Teak (Tectona grandis Linn. F.) is one of the most valuable tropical hardwoods of the world, being sought in the global markets for its beauty, strength and stability, natural resistance and wide array of applications [1]. The importance of teak plantations around the world can be appraised from the report "State of the World's Forest Genetic Resources" [2]. It lists tree species that are considered national priorities by the reporting countries for the conservation and management of forest genetic resources, teak taking the top rank in this list in more than 20 countries due to its economic value.
In 2010 the global area of planted teak forests reported from 38 countries was conservatively estimated at 4.35 million ha, which 83% grew in Asia, 11% in Africa, 6% in tropical America, and less than 1% in Oceania [3]. More recent estimations of the area of planted teak forests based on FAO (Food and Agriculture Organization of the United Nations) and ITTO (International Tropical Timber Organization) reported up to 6.89 million ha in 2015 [4].
Regarding Ecuador, and according to official data provided by MAGAP (Ecuadorian Ministry of Agriculture, Livestock, Aquaculture, and Fisheries), it currently has around 50,000 ha of planted teak forests, mainly located in the so-called Coastal Region (provinces of Guayas, Manabí, Esmeraldas, El Oro, Los Ríos, Santo Domingo, and Santa Elena). Note that MAGAP fostered the establishment of pure teak plantations in coastal and Amazonian Ecuador within a governmental program of incentives for reforestation with commercial purposes launched in 2013, establishing a goal of 120,000 ha of new planted teak forests in five years. However, progress in the suggested goals expressed as the reforested area has been modest [5]. In any case, the economic value of teak plantations for the country of Ecuador is supported by the figures of average global annual teak roundwood export volumes between 2005 and 2014, which amounted to values between 50,000 and 100,000 m 3 /year [6].
Forests are, in general, net carbon sinks, thus making forest conservation of critical importance for mitigating climate change [7]. Since forests can conserve or sequester significant quantities of carbon [8], slowing deforestation, combined with an increase in forestation and other management measures to improve forest ecosystem productivity, could foster negative emissions and contribute to meet the Paris agreement aspiration of limiting global temperature rises to below 2 °C [9]. This sequestration of atmospheric carbon by forests also becomes a major strategy for the United Nations Framework Convention on Climate Change within the context of Reducing Emissions from Deforestation and forest Degradation (REDD), which could help mitigate greenhouse gas emissions on forest-rich developing countries [10], as is the case of Ecuador. In fact, it is widely recognized that reducing carbon emissions from deforestation and degradation in developing countries is of central importance in efforts to combat climate change, which requires information on forest clearing and carbon storage [11].
Despite the need for increasingly precise and efficient forest monitoring, characterization of forest plantations at tree level has been limited to field-based techniques focused on measurement of diameter of the cylindrical part of the bole and other size attributes such as total tree height and crown size, involving large uncertainty in measuring large trees with irregular shapes [12]. In this way, traditional studies aimed at developing allometric models to estimate dry above-ground biomass (AGB), or other tree-level variables, are based on costly and inefficient destructive sampling methods, trying to construct allometric models applicable worldwide [13,14] which generally lead to uncertainties and systematic errors in biomass estimates when applied under local conditions.
In response to the increasing demand for efficient tools for forest monitoring, a new paradigm has emerged that cannot be seen as a logical progression of the measurement based on existing plots, but as a completely new way of dealing with forest inventories that highlights the use of remote sensing technologies [14][15][16][17]. That is the case of terrestrial laser scanning (TLS), an efficient and nondestructive measurement method that should be viewed as a disruptive technology that requires a rethink of vegetation surveys [18][19][20][21][22][23][24] to support the development of locally adjusted allometric models [25][26][27][28] and the 3D description of tree architecture [29][30][31][32]. This detailed 3D tree description can be used for the estimation of AGB and carbon stocks via volume estimation through using the so-called quantitative structure models (QSMs) [27,33].
The main objective of this study aims at developing locally calibrated allometric equations based on TLS data obtained in leaf-off conditions teak plantations located in the Coastal Region of Ecuador. Together with the unmanned aerial vehicle structure from the motion (UAV-SfM) method tested in a previous work [34], it could make up the skeleton of an efficient remote sensing based method to estimate and update the aboveground biomass and tree stem commercial volume dynamic over small sample areas of the teak plantations located in the study area (i.e., at project level scale within the context of REDD), thus contributing to their sustainable management and carbon stock estimation.
The study site presents an average rainfall of about 1222 mm, with an average temperature of 24.4 °C and a relative humidity of 72.9% [35]. It belongs to the Tropical Dry Forest [36], being characterized by an unimodal rainfall regime with a rainy period in the first quarter of the year and a marked drought during the rest of the year. These climate conditions cause alternating leaf-on (April to September) and leaf-off (November to December) seasons in the case of the teak plantations located at the coastal region of Ecuador. The data on which this research is based were collected during the leaf-off canopy conditions in order to facilitate the field work avoiding shadows and occlusions [34].

Field Data
The forest inventory based on processing terrestrial laser scanning (TLS) point clouds offers, in addition to a very precise position of each tree (absolute georeferencing) and the extraction of dendrometric variables at tree level such as diameter at breast height (DBH) and total tree height (h), very valuable geometric information about the complete geometry of each tree stem, which allows even the estimation of its stem taper and quasi-real volume [19][20][21]37].
A TLS field campaign was carried out in November 2018 (leaf-off conditions) over 58 reference plots of even-aged teak trees according to the following distribution: Thirty plots were located in Morondava teak plantation, eight in El Tecal plantation, and the remaining 20 in Allteak plantation (see location in Figure 1). Some key features of the reference plots located in the plantations of Morondava, El Tecal, and Allteak can be found in Aguilar et al. [34].
Morondava would be a good representative of a typical young teak plantation with plantation ages between two and three years (2.9 years on average) and an average planting density (surviving trees) of 722 trees/ha. This is a very early developing plantation that shows a certain heterogeneity with both a low average basal area of 4 m 2 /ha (from 0.6 to 6.28 m 2 /ha) and a low average Lorey's mean height of 7.83 m (from 3.84 to 9.57 m). The average slope of the reference plots was 23.4%, ranging from 3.9% to 55.5% [34]. Understory vegetation was mainly composed of shrubs (Mimosa pudica Linn.) and an herbaceous layer (Convolvulus arvensis Linn.). However, ground level was detectable since herbaceous cover was relatively transparent at the time of field work. Clearing was needed to remove shrubs and some herbaceous cover before carrying out the TLS field work.
Regarding El Tecal, it is a 17-year-old plantation much more developed than Morondava and with an average planting density of 985 trees/ha. Both basal area and Lorey's mean height presented average values of 16.78 m 2 /ha and 15.36 m, respectively, significantly higher than those registered in Morondava. Specifically, the average basal area increment with respect to Morondava was 319.5%. Considering that both basal area [38] and Lorey's mean height [39] are dasometric variables positively correlated to forest stand volume, it can be inferred that the reference plots located in El Tecal presented a forest stand volume larger than those of Morondava. The average slope of the reference plots was 11.7%, ranging from 8.3% to 14.7% [34]. Understory vegetation was mainly composed of an herbaceous layer (Amaranthus retroflexus L. and different species of grasses). Shrubs were scarce due to high competition with teak trees. Ground level was practically covered with fallen teak leaves.
In the case of the reference plots located in the Allteak plantation, the age of plantation varied between four and 12 years, showing a very heterogeneous planting density with an average of 579 trees/ha, a value lower than that registered in Morondava or El Tecal due to the usual application of thinning operations. Basal area values ranged from 5.17 to 17.12 m 2 /ha (11.56 m 2 /ha on average, meaning an average basal area increment with respect to Morondava of 189%), while Lorey's mean height values varied from 12.30 to 23.02 m (18.57 m on average). The average slope of the reference plots was 16.8%, ranging from 2.5% to 43.9% [34]. Understory vegetation was composed of a dense herbaceous layer (grasses such as Rottboellia exaltata L. and Panicum maximum L.) and a scarcely developed shrub layer due to competition with teak trees. Ground level was entirely covered with fallen teak leaves. Clearing was needed to remove herbaceous cover before TLS field work.
A very dense and accurate TLS point cloud was obtained within each reference plot through a field survey carried out with a FARO Focus 3D X-330 TLS instrument (FARO Technologies Inc., Lake Mary, FL, USA), collecting x, y, and z coordinates together with high resolution RGB images which served as a high quality visual reference in post-processing tasks. The TLS device, which utilizes phase-shift technology, was set to acquire data with medium resolution and quality (1/5 resolution and 4× quality) for potentially obtaining up to 28.2 million pulses per scan along a range from 0.6 to 330 m and a vertical and horizontal field of view of 300 and 360°, respectively. It is a portable TLS device (weight of 5.2 kg) that works with a wavelength of 1550 nm and a beam divergence of 0.19 mrad, yielding a ranging error of ±2 mm (10 to 25 m).
Four scanning positions were set up within a radius of 18 m from the center of each reference plot to configure a scan pattern with a central scan and the rest located around, drawing an equilateral triangle ( Figure 2a). It has been demonstrated that multi-scan data (merged scans) are more accurate to extract stem diameter and stem volume than single scans [40]. The FARO Scene © 7.1 software (FARO Technologies Inc., Lake Mary, FL, USA) was utilized to co-register the four scans within each reference plot, thus producing a single and colored 3D point cloud, from using nine artificial targets (15 cm diameter spheres) conveniently distributed over the reference plot to ensure that at least three spheres were visible from every two consecutive scan positions. The mean distance error of the co-registered point clouds for the 58 reference plots was 6.2 mm, with maximum and minimum values of 13.3 and 2 mm, respectively.
The number of points per square meter was variable depending on the vegetation coverage of each reference plot. In the case of Morondava plantation, the point density ranged from 19,771 to 34,776 points/m 2 (28,472 on average). Point density ranged from 28,090 to 33,868 points/m 2 (31,859 on average) in El Tecal plantation, while it varied between 29,472 and 63,038 points/m 2 (40,820 on average) in Allteak plantation.
The scan located at the center of each reference plot (reference scan for the co-registering process) was georeferenced by applying a 3D conformal coordinate transformation based on the WGS84 UTM 17S coordinates of four spheres (out of the nine available) surveyed by using a TOPCON ES105 total station. It was oriented for each reference plot from the observation of up to three GPS-RTK points, located very close to each reference plot, that were previously surveyed during a field campaign conducted between March and June 2018. Time elapsed per scan was approximately 8 min, while completing a reference plot took around 93 min on average (including clearing of understory vegetation, spheres positioning, and absolute georeferencing).
A very accurate scattered bare earth points were obtained from automatically segmenting TLS ground points by means of the octree search algorithm implemented in the open-source software 3D Forest [41] (Figure 2b). These bare earth points allowed to build a highly accurate 20 cm grid spacing digital terrain model (DTM) as ground reference to compute the heights (normalized heights) of every point belonging to the TLS point cloud.

Tree Database Collection
The pipeline shown in Figure 3 was applied to the plot-level collected TLS point clouds in order to attain a representative database of teak trees from which extracting dendrometric variables such as DBH, h, tree stem volume, and tree stem commercial volume.

Automated Tree Segmentation
The open-source software 3D Forest version 0.42 (www.3dforest.eu) (accessed on 9 September 2019) was used to cope with the segmentation of forest vegetation (i.e., previously classified nonground points) into individual trees. Segmentation is based on point clusters, which are defined by maximal distance between the two nearest points in the cluster (S (cm)), minimal number of points forming clusters (N), and the angle (10°) and distance between centroids (4S) of the clusters to add new clusters to the previously reconstructed bases of the trees. The reader can find a detailed description of this algorithm in [41]. The corresponding tuning parameters were tested before obtaining the automatic segmentation result in every reference plot. The best value for parameter S was 15 cm in all plantations and reference plots, while the optimal value of N varied from 5 to 30 points (greater when the slope of the reference plot was greater). Finally, an additional manual editing was needed to correct segmentation errors by deleting redundant parts from the wrong tree and adding them to the right one. Note that this type of automatic segmentation algorithms perform better in the leaf-off phenological stage due to the presence of fewer obstacles (thus fewer shadows in the point cloud) [41,42]. In this sense, we fully recommend carrying out the TLS field work in leafoff conditions (dry season in the case of the Coastal Region of Ecuador).

Selection of Suitable Trees
Up to 3791 teak trees (i.e., TLS segmented point clouds) were extracted from the 58 reference plots previously described. After a data filtering process to remove trees with DBH < 5 cm and/or crown diameter less than 1 m, which corresponded to small trees that showed poor crown development at the time of scanning (mostly belonging to Morondava teak plantation), 2272 trees were chosen to carry out the robust regression model that relates DBH (dependent variable) and h (predictive variable).
From the total of 2272 teak trees that made up the original database, 456 were selected because they presented a suitable geometry according to criteria such as no bifurcated trees, low inclination, or absence of epiphytic vegetation. In addition, the selected trees were intended to cover the widest possible range of the biophysical variables of the teak plantations where the study was conducted, trying to count on a balanced distribution of sampled trees according to their DBH (Table 1). The 456 selected teak trees were manually edited using the free and open source software 3D Forest (https://www.3dforest.eu/#download; accessed on 9 September 2019) to filter out points that did not belong to the tree stem such as leaves and branches, allowing the isolation of the TLS points that formed the tree stem geometry ( Figure 4). Note that although there are algorithms capable of isolating the stem from the rest of the parts of the tree (e.g., [43]), the goal of this work was not based on testing or improving this type of algorithms but on obtaining reliable allometric models. This is the reason why a careful and time-consuming manual edition was carried out in order to ensure that the tree stem (from the stump to the top) was as geometrically clean as possible.  The same segmented tree after manual editing from which the points have been manually classified in stem tree points (red color) and branches and leaves tree points (green color).

Automatic Extraction of Tree Level Information
A MATLAB code software called Tree_geometry was developed for the automatic extraction of the variables of interest at tree level such as DBH, h, and commercial and total stem volume.
Tree_geometry makes cross sections along the tree stem at different heights (hi) (normalized heights with respect to the ground), starting with two sections at heights of 0.15 and 0.30 m. From the height of 0.30 m, sections are produced according to increments of 0.50 m, thus including in the dataset the value of 1.30 m where the DBH is measured. Each cross section hi encompasses points from the tree stem TLS point cloud located at a height hi ± 2.5 cm. The section is considered valid if the number of points is greater than 20. This means that there may be some missing cross sections because they have few or no TLS points due to occlusion problems during the scanning process. The points included in each cross section i, after their 2D projection onto the horizontal plane, are adjusted to a circle using the robust adjustment method proposed by Ladrón de Guevara et al. [44]. This method of adjustment is very suitable when outliers are foreseen due to TLS points not belonging to the tree stem (e.g., those located in epiphytic vegetation, nearby branches, or understory vegetation). It is also recommended to deal with the adjustment of incomplete circles (circle arcs), as can be the case in TLS inventories when the entire circle of a tree stem cross section is not covered because of occlusion problems or simply due to an insufficient number of scanning positions in the reference plot [19]. Since most geometric or algebraic adjustment methods are sensitive to noise and atypical points, Tree_geometry uses the minimum absolute error (MAE) criterion to provide robust estimates. Considering that the MAE objective function is not differentiable, and so the adjustment methods based on the gradient cannot be applied, the implemented algorithm relies on the partial derivative on the right and on the left (exterior and interior points to the reference circle in each iteration). It allows to significantly increase its computational efficiency, representing a good alternative to those methods based on least squares fitting, RANSAC-based least squares fitting, or randomized Hough transform [18,22,41].
The commercial stem volume has been computed from a height of 0.15 m, which turns out to be the usual stump height at which the teak tree is usually cut [45]. The commercial height was determined as the height of the tree stem with a diameter of 0.13 m (measured with bark), i.e., the minimum usable diameter according to the usual commercialization criteria [46,47]. The cases in which the first bifurcation occurred before reaching a diameter of 0.13 m were filtered out through the selection procedure of the teak trees that make up the database (see Section 2.3.2). The commercial height was linearly interpolated between the height values of the 0.5 m separated cross sections whose corresponding diameters were immediately above and below the target diameter of 0.13 m. The volume of each log (volume Vi,i+1 between two consecutive cross sections Si and Si+1 of radii Ri and Ri+1) was computed as the volume of a truncated cone given by the following expression: where hsi,i+1 being the difference in height between the two consecutive cross sections Si and Si+1 (usually 0.5 m). The sum of the volumes of all the logs, including the final one up to the diameter of 0.13 m, would be a good estimation of the commercial stem volume of each tree. The total tree stem volume was determined by adding two components to the commercial stem volume: i) The stump volume, calculated as a cylinder of 0.15 m height and the diameter corresponding to the cross section located at 0.15 m height; ii) the volume of the final apex of the tree stem located between the highest cross section of the every 0.5 m sections and the maximum height of the tree stem. This final apex was modeled as a cone.

Development and Validation of Allometric Models
Several mathematical formulations of allometric models based on the widely known potential form (e.g., [14]) were tested. They used the following general expression: being Y the dependent variable to estimate (commercial and total stem volumes in our case), X the independent or predictor variable at tree level (DBH, tree height, and combination of both), and a and b regression model coefficients. The parameters of the allometric model were determined by applying regression method after a previous logarithmic transformation of the variables for adjustment using a linear model. In this way, the following allometric models were explored [13,14,48]: lnY =∝ +β ln(h) + ε, where DBH is expressed in centimeters and h in meters. α and β are model coefficients (derived from least squares regression), and ε is an error term, which we assume to be normally distributed with zero mean and standard deviation σ [13,14]. It is relevant to highlight that the formulation of Equation 3 turns out to be the well-known allometric model proposed by Pérez and Kanninen [49] to be applied in teak plantations. This allometric model estimates the total AGB (kg/tree) according to the following expression: Considering the bias that generally occurs through the least squares estimation of the regression parameters, robust regression is recommended. This bias is motivated by the unbalanced distribution of tree size in the dataset (small stems are generally more frequent) and the likely presence of outliers [14]. In this way, a robust regression procedure has been programmed using MATLAB code, developing an iterative algorithm for assigning weights to data using a bi-square weights function [50,51]. The bi-square weights function adopts the form indicated in Equation 7, where ri represents the value of the normalized residue for the i-th point which, in turn, is given by the expressions shown in Equation 8. The validation of the tested allometric models was based on the so-called true validation method, a validation procedure meaning that the data used to calibrate the model were never used for its validation. In this sense, the validation set consisted of 10% of the available trees (from the 456 trees available. See Table 1), leaving the remaining 90% as a set for training and computing the regression coefficients. This procedure was repeated up to 100 times, performing a data extraction and assignment to the validation set employing a stratified random sampling based on the distribution of values of the DBH. In order to validate the results and obtain the uncertainty of the regression model, the RMSE (root mean square error), relative RMSE (estimation accuracy), Bias (systematic error), and σ (residual standard error or regression uncertainty) values were computed on every validation dataset extracted by stratified random sampling according to Equations 9, 10, 11, and 12, respectively [14]. The robustness and variability of the estimates obtained for the RMSE, Bias and σ values could be determined from the inherent variability of the 100 validation datasets available.
where n is the number of observations in the validation dataset, and Ypred and Yobs refers to the predicted and observed values, respectively. The final allometric models can be used to predict the dependent variable (commercial or total stem volume) supposing that DBH and/or h are known. For example, the predictive model for the allometric model based on only DBH (Equation 3) would be as follows: Ypred = e (α+βln(DBH)+ε) . As it was said before, ε is assumed to be normally distributed [i.e., N(0, σ 2 )], therefore the mean of can be approximated by [52]. Note that the term can be understood as a correction factor applied to back-transform predicted values and remove bias from the log-transformed data [28,52]. An unbiased estimate of Ypred can therefore be calculated using the following equation: Since the relationship between DBH and h is the most commonly used measurement of tree size [53], the knowledge of this relationship is especially relevant to develop allometric models mainly based on tree attributes, such as tree height, which can be remotely sensed from spaceborne and airborne sensors in order to increase the extension of the remotely sensed forest area [14,16,34,54]. In this way, a power-law relationship between h and DBH (e.g., [55]) was tested to estimate DBH from total tree height observed values according to a mathematical model such as that shown in Equation 13 (specifically ). The inclusion of the crown diameter parameter as an additional explanatory variable was ruled out due to difficulties in an accurate measurement of the crown size from the TLS data collected in leaf-off conditions.
The sample of 2272 teak trees previously described in Section 2.3.2 was used to fit the model DBH = f(h) by applying the described robust regression method according to two different strategies: • Global fitting. The model parameters were computed from all available sample teak trees (2272 individuals).
In addition, and due to the likely presence of outliers in the measurements taken from the original dataset, some additional robust statistics such as median (M) and median absolute deviation (MAD) were computed over the residuals given by the difference between the predicted and observed values for both commercial and total tree stem volumes. The median is a measure of central tendency that offers the advantage of being very insensitive to the presence of outliers. The same can be said about MAD as an estimator of scale or variability. Moreover, MAD is not sensitive to the sample size [56]. MAD was calculated from the following expression: where xi refers to the set of n observations (residuals in our case), M is the median of the observations and k is a constant scale factor, which depends on the distribution. Usually k = 1.4826 if normality of the data is assumed, disregarding the abnormality induced by outliers [57].

Relationship Between DBH and Total Tree Height
The validation results of the global fitting that relates DBH and h are depicted in Table 2 (column 2). This model exhibited a relatively low RMSE (1.69 cm; 13.85% relative RMSE), slightly overestimating the observed DBH values with a bias of 1.32%.
On the other hand, the DBH-based fitting performed better than the globally fitted model ( Table  2, columns 3 to 5), providing lower RMSE values, especially for the greater range of DBH. This improvement was motivated because the DBH-h relationship turned out to be nonlinear on a log-log scale for the entire DBH range, as described earlier in some studies [13,14]. from all available sample trees (Global fitting) and from sample trees grouped by DBH classes (DBH-based fitting). The standard deviation calculated on the 100 validation datasets is given in brackets. Since the main objective was to apply the DBH-h relationship for predicting tree stem volume and AGB values at tree level, we tested replacing DBH with both the locally and globally fitted models in order to check the corresponding error propagation through the allometric model proposed by Pérez and Kanninen [49] (see Equation (6)).

Statistics
The results can be seen in Table 3, depicting how the globally adjusted model produced highly biased AGB predictions. Indeed, it overestimated the AGB for small trees (that is, low DBH values) presenting a bias of 23.24%, while clearly underestimated the AGB in the case of trees with DBH between 20 and 30 cm (bias of -19.84%). In addition, the relative RMSE reached values of up to 47.63% in the case of the estimation of the AGB of small trees, noting that the globally fitted DBH-h model was not adequate to be incorporated into the development of allometric equations to predict AGB. On the contrary, the DBH-based fitting model was able to produce relatively unbiased estimates of AGB across the whole range of DBH, also keeping the relative RMSE values below 28%. This was caused by a smoothing of the extreme AGB predicted values, both the lowest and the highest, due to a better adjustment of the log-log DBH-h relationship in the higher and lower DBH ranges. It is worth noting that this precision gain is of high interest to reduce the bias of allometric models when DBH is replaced by h.

Allometric Model to Estimate Tree Commercial Volume
After applying the robust regression method previously described, the three allometric models shown in Table 4 were obtained. Furthermore, another allometric model only based on tree height was derived from replacing DBH by h in the expression DBH 2 × h (Equation (5)), according to the DBHbased fitting model previously described.  It can be seen how all the primary models presented a positive bias between 8.64% to 12.72%, which indicates an overestimation of the observed commercial stem tree volume. The RMSE had values close to 0.05 m 3 , which implied an estimation error (relative RMSE) of between 21.6% and 33.22%. The two-variable model, which had DBH and h as explanatory variables, turned out to be the most accurate, with relative RMSE and RMSE values of 21.6% and 0.0482 m 3 , respectively. The allometric model based only on tree height, an attractive model because tree height is relatively easy to obtain by means of UAV-based stereo-photogrammetry or airborne lidar [12,58], presented an estimation error of 33.22% (RMSE = 0.074 m 3 ), so it would allow a reasonable approximation for the estimation of the commercial volume of teak plantations.
However, the high bias of the allometric model based only on tree height (12.72%) recommended testing the results provided by replacing DBH with h in the two-variable allometric model, which was the primary model that offered the best validation results. The results of this derived model (Table 4) notably improved those obtained with the model based only on tree height, decreasing its bias to only 1.43%. Furthermore, its relative RMSE was 27.24%, while its RMSE turned out to be 0.0608 m 3 .
The goodness of fit of the derived allometric model only based on tree height can be graphically observed in Figure 5. It is worth mentioning the presence of possible outliers that have been adequately excluded from the adjustment process thanks to the robust regression method applied, together with a marked evidence of underpredicting when the observed commercial volume values were above 0.5 m 3 . Despite these shortcomings, when robust statistics insensitive to outliers such as median (central tendency) and MAD (variability) were computed over the residuals, reasonably low values of −0.0103 and 0.0407 m 3 were obtained, respectively. It should be highlighted that, in this case, the allometric model only depends on tree height, a variable that can be remotely sensed from spaceborne or airborne sensors.

Allometric Model to Estimate Tree Stem Volume and Dry Biomass
The results of the allometric models tested to predict tree stem volume are shown in Table 5. The two-variable model (DBH 2 × h) is highlighted since it performed predictions with a practically negligible bias (overestimation of 0.37%), a low RMSE of 0.0319 m 3 , and a relative estimation error of only 16.41%. The allometric model based on only tree height performed slightly worse, presenting a bias of 2.54% and RMSE and relative estimation error of 0.0503 m 3 and 25.88%, respectively.
Although the results provided by the allometric model based on only tree height were acceptable, and as in the previous section, a derived model was tested. As said before, this derived model consisted of replacing DBH in the two-variable allometric model by using the DBH-based fitting model developed in Section 3.1 (Table 5, fourth model). The results provided by this derived model were slightly better than those of the allometric model based only on h, presenting a bias of -1.11% and a RMSE of 0.0461 m 3 , which corresponded with a relative estimation error of 23.72%. In this sense, this derived model would allow predicting tree stem volume with an average error of less than 25% in teak plantations only from knowing tree height. Figure 6a shows the plot of observed/predicted tree stem volumes for observed values ranging from very low values to around 0.65 m 3 . As in the case of commercial volume, there was a tendency to underpredicting tree stem volume at the largest tree sizes, that is tree stem volumes above 0.5 m 3 . In addition, it is worth noting the likely presence of outliers in the dataset since very low values of −0.002 and 0.0247 m 3 were obtained when computing robust statistics such as median and MAD, respectively. Table 5. Allometric models to estimate tree stem volume (Vstem in m 3 ). The explanatory variables are DBH (cm) and/or h (m). To the best of the authors' knowledge, currently there is only one allometric model developed in the Coastal Region of Ecuador. This is Lara's model [59], which was locally fitted to estimate tree stem volume in teak plantations located in Balzar (province of Guayas) using a traditional destructive sampling. The sample was composed of 125 teak trees with DBH between 1.8 and 40.8 cm and tree heights between 2.1 and 22.7 m. Lara's model depends on both DBH and tree height to compute the tree stem volume according to the following expression: where Vstem is given in m 3 , and DBH and h are expressed in cm and m, respectively. The fit of this model to the sample data collected in this work is shown in Figure 6b, pointing to a clear underestimation of the actual values even for medium size trees. In quantitative terms, the bias of the adjustment was −6.80%, presenting a RMSE value of 0.0521 m 3 and a relative RMSE of 26.79%, so performing worse than any of the allometric models developed in this investigation (see Table 5). This performance difference was especially notable in the case of the allometric model based on the metric DBH 2 × h, a comparable model because it uses the same input variables, that provided unbiased results (bias = 0.37%) and low RMSE and relative RMSE values of 0.0319 m 3 and 16.41%, respectively. .
,where f(h) represents the DBH-based fitting model which relates DBH and h. b) Allometric model proposed by Lara [59]. The red line refers to the 1:1 line.
The tree stem dry biomass was estimated through converting the tree stem volume in dry biomass by applying a dry density value of teak wood of 550 kg/m 3 , an average value extracted from bibliography for areas such as Ecuador, Brazil, Venezuela, Bolivia, and Costa Rica [60][61][62][63][64]. It should be borne in mind that the density of teak wood is relatively variable with the age of the tree, increasing with increasing age [62].

Allometric Model to Estimate the Above-Ground Biomass at Tree Level
In order to estimate the dry above-ground biomass (AGB) at tree level, it is necessary to add the contribution of branches and leaves to the biomass previously estimated for the tree stem. In this sense, a fixed ratio between the total AGB and the biomass contributed by the tree stem was adopted. The results published by Pérez and Kanninen [49], obtained in teak plantations in Costa Rica through the destructive sampling of 87 trees with ages between five and 47 years and DBH range between 10 and 59 cm, reported variations in the relative distribution of dry tree stem biomass between 70% and 90% with respect to the total tree AGB. According to these results, an average value of 80% was applied in this work, so that the total tree AGB was estimated as the dry tree stem biomass divided by 0.8. Note that this is only a mere approximation that should be adapted in each case, being beyond the scope of this work to reconstruct and modeling complete trees from TLS point clouds, a very active and promising research line still under development (e.g., [29,32]). Our approach will allow the user to estimate the total AGB of each teak tree based on the apparent crown size volume, so rating its contribution to the total biomass and being able to modify, if deemed necessary, the factor of 0.8 used in this work.
In this way, the allometric model of AGB at tree level would be given by the following expressions: • If we only have the tree height: Equation 15 would be applied, but substituting DBH for its estimation from the following expressions (see the DBH-based fitting model parameters in Table 2 (20 cm ≤ DBH < 30 cm) Figure 7 shows a graphical comparison between the allometric models for the estimation of AGB proposed by Pérez and Kanninen [49] (Equation 3), based only on DBH measurement, and the two AGB allometric models developed in this work (Equations 16 and 17). The great similarity between the results provided by the two-variables (DBH and tree total height) allometric model as compared to the Pérez and Kanninen AGB allometric model can be seen in Figure 7a. If the model of Pérez and Kanninen is taken as a reference, the two-variables based model performed quite similar, yielding unbiased AGB values (bias = −0.12%) with low RMSE and relative RMSE values of only 18.27 kg/tree and 14.41%, respectively. The goodness of fit between both models proves the soundness of the TLS-based methodology (non-destructive method) proposed in this work to efficiently build locally adapted allometric models in teak plantations.
When we tried to avoid the need of DBH, also taking the Pérez and Kanninen model as reference, the comparison results were clearly worse (Figure 7b). In fact, the AGB allometric model based only on tree height, despite being very useful due to the ease of measuring tree height from airborne or spaceborne sensors, yielded slightly biased results with an AGB underestimation of −1.38% with respect to Pérez and Kanninen model. However, the proposed model clearly underestimated the reference of AGB values in the case of larger trees (more than approximately 350 kg/tree), showing a median of AGB differences between the two models of −3.49 kg/tree. In the same way, the variability of residuals was higher than in the case of the two-variables AGB model, obtaining an RMSE value of 34.17 kg/tree, which meant a relative RMSE of 26.95%, and a MAD value of 19.70 kg/tree.  Table 6 depicts the average and standard deviation values of each DBH class in which we grouped the sample of teak trees with respect to dendrometric variables such as tree stem commercial volume, tree stem dry biomass, and AGB. It should be noted that the tree stem commercial volumes obtained in this work coincided with those reported by Krishnapillay [65], who provided values of 0.5 m 3 on average. Similarly, Mora and Hernández [66], working on teak plantations in the Pacific of Costa Rica, reported an average commercial volume of 0.48 m 3 for the DBH class between 25 and 30 cm, a value that fully coincides with that obtained in the present investigation. Table 6. Average values of some dendrometric variables grouped by DBH classes. The standard deviation of the observations is noted in brackets. Sample size = 456 trees (see Table 1).  Table 6 are only intended as a guide and reference for further researches.

Discussion
It is necessary to highlight the added value of the methodology and allometric models provided by this study for the sustainable management of teak plantations located in the Coastal Region of Ecuador. In fact, with the exception of the allometric model proposed by Lara [59], there are no local allometric models available in this area to estimate AGB at tree level, being the model published by Pérez and Kanninen [49], which was developed in teak plantations located in Costa Rica, commonly used (see Equation 6).
On the contrary that Lara's allometric model, which requires both DBH and tree height as inputs, one of the strengths of the allometric models proposed in this work is that they can be applied only knowing the tree height from the mathematical relationship found between tree height and DBH (Table 2), especially the one based on DBH ranges (DBH-based fitting) that showed a reasonably low error propagation to calculate AGB from applying the allometric model of Pérez and Kanninen ( Table  3). Note that current allometries rely on DBH as a key input for estimating AGB [13], but stem diameter cannot be measured directly from airborne or spaceborne sensors, which are the most suitable remote sensor based technologies to carry out wall-to-wall forest inventories by determining tree height from the widely known canopy height model (CHM) [16]. Counting on allometric models based solely on tree height could allow to compute the average tree height within the area of interest using the CHM [34,67,68], while the number of trees per hectare could be approached simply applying a maximum local filter over a variable window size [69] or more sophisticated methods based on point cloud 3D clustering [12]. This would allow to conduct an agile area-based approach (ABA, a distribution-based technique which typically provides data at the stand level) to estimate AGB or even commercial volume in teak plantations.
The tree height based allometric models developed in this study allowed to estimate the tree stem commercial volume, the total tree stem volume, and the tree stem biomass, showing a low bias and a reasonably low relative RMSE. However, the underestimation of the values observed in the case of larger trees is remarkable, probably due to a failure of the power law ratio h-DBH to adequately capture the asymptotic nature of tree height growth [13,70,71]. Trees, especially fast-growing trees such as teak, generally invest primarily in height growth when they are young, rapidly approaching their maximum height, but then continue to grow in diameter throughout their lives [14,46]. This makes linking the DBH with the tree height challenging in the case of large trees, since trees of similar height can have very different DBH values. A complementary explanation of the underestimation of tree stem volume and AGB predicted values in larger trees could be related with the underestimation of tree height. In fact, the main source of underestimation of tree height by TLS is that the top of the tree crown is occluded by itself or a neighboring tree [18]. This drawback was partially resolved by carrying out the field work in leaf-off conditions, that is, during the dry season of the Coastal Region of Ecuador. However, some plots, especially those located on the Allteak plantation, which housed the largest trees in the sample, still contained trees that generally had dry leaves that had not yet fallen to the ground since they were in areas with greater availability of soil water due to their edaphological and/or geomorphological characteristics. For example, they were located in depressions or channels where rainwater accumulates. This circumstance probably caused an underestimation of tree height which, in turn, affected the adjustment of the allometric models in the case of the largest trees in the sample.
It is widely known that while height growth tends to slow down rapidly in large trees, lateral crown expansion does not, requiring a continued investment in stem growth to ensure structural stability and hydraulic function [14,71]. Consequently, crown width and stem diameter tend to be strongly coupled, even in large trees [72]. In this regards, we tried to include the average crown diameter (CD) in our model to estimate DBH by using the metric h x CD by means of the mathematical model = (∝ ( . )) [14], but the results were not significantly better than those obtained from the model only based on tree height (data not presented in this work). It was probably due to an insufficiently accurate measurement of CD from TLS data in the case of leaf-off conditions, together with the distortion effect of pruning on the natural relationship between DBH and CD.
As explained in Section 3.4, a fixed ratio of 0.8 between the biomass contributed by the tree stem and the total AGB at tree level was adopted, considering the results published by Pérez and Kanninen [49]. This approach turns out to be very simple, being necessary to improve the estimates of the biomass contributed by the branches from TLS-derived point clouds (non-destructive method) as a future research direction. In fact, reducing uncertainties in tree and forest carbon estimates requires an accurate AGB prediction [11]. Accordingly, two alternative methods could be faced.
The first one relies on the promising and emerging research line related to the automatic or semisupervised modelling of laser-scanned trees including tree stem and branches (e.g., [29]). Note that leaves are not usually included in tree modelling, since it is almost impossible to capture the geometry of leaves from laser scans, therefore, becoming impossible to reconstruct accurate leaves purely from the point clouds [20].
The second method is based on tree physiology, considering that continual new-leaf production is essential for the maintenance of photosynthesis. In response to new-leaf production (primary growth), stem diameter growth (secondary growth) would occur mainly within the crown. This would allow newly produced leaves to secure conductive tissues [73]. It is therefore likely that the priority of growth within the crown corresponds to the need for leaf turnover [53], supporting the hypothesis that the stem diameter measured at a certain height within the crown could be a good indicator of the quantity of branches and leaves of a tree. It is worth noting that this hypothesis would coincide with the pipe model theory which proposes that the square of the diameter of the stem just below the crown base (crown-base diameter) is directly proportional to the amount of leaves on an individual tree [74]. In this sense, 3D modelling of teak trees from TLS-derived point clouds would allow to find local relationships between crown-base diameter (or whatever stem diameter within the crown) and the volume of tree branches.

Conclusions
The results obtained in this study have proved that TLS-derived point clouds constitute an efficient tool to capture teak tree geometry at stem level under leaf-off conditions in the Coastal Region of Ecuador, although further work is needed to deal with an accurate extraction of tree branches in order to count on the full tree skeleton.
The allometric models, calibrated and validated for planted teak forests located in the Coastal Region of Ecuador, provided estimates of some tree-level variables, such as tree stem commercial volume, tree stem volume, and AGB, with a reasonably low relative error (relative RMSE) and almost no bias. In addition, it is necessary to highlight the development of a mathematical relationship between DBH and total tree height that allows the estimation of DBH from total tree height based on expressions which were adjusted according to the range of DBH (specifically three classes of DBH). This set of tools can facilitate the obtaining of variables such as AGB (sequestration of atmospheric carbon by teak plantations) and commercial volume of teak wood (economic variable) in the Coastal Region of Ecuador from the knowledge of simply the total tree height. This finding constitutes a significant advantage to address digital inventories based on airborne sensors such as unmanned aerial vehicle-based structure from motion (UAV-SfM) method, enabling large-scale forest inventories (upscaling) through production of accurate canopy height models.
The allometric equations developed in this study, together with the UAV-SfM method tested in a previous work [34], make up the basis of an efficient remote sensing based method to estimate treelevel forestry parameters, such as AGB or tree stem commercial volume, over teak plantations located in the Coastal Region of Ecuador. This method could potentially be applied, after proper calibration and validation, to other teak plantations located in seasonally dry tropical forests to assist in their REDD monitoring and sustainable management.
Author Contributions: F.J.A. provided the overall conception of this research, designed the study, and developed the software, also writing the original manuscript and participating in the data processing. A.N. conducted the TLS fieldwork and processed TLS data, also helping in the experimental design and general data processing and curation. A.P. conceptualized and supervised the sampling design and the forest inventory fieldwork, also contributing to the project administration, funding acquisition, and overall conception of this research. J.R.R. conducted the field work, also helping in providing the resources needed to carry out the forest inventory. M.A.A. contributed to outline the experimental design, statistical analysis, and manuscript conception. A.P., A.N., M.A.A., and J.R.R. revised the original manuscript. All the authors analyzed and discussed the experimental results together.
Funding: This research was funded by the project "Evaluación de tecnologías de detección remota para la estimación de biomasa de Teca en la Región Costa de Ecuador" (Research and Development System of the Catholic University of Santiago de Guayaquil, Ecuador). This work also takes part of the general research lines promoted by the Agrifood Campus of International Excellence ceiA3, Spain (http://www.ceia3.es/en).