Introduction

The state of Quintana Roo (Q. Roo) occupies the second national place in timber production of broad-leaved species, with an approximate harvest of 47 139 m³ of roundwood (m³ mr) of common tropical species (12.10 % of national production) and of 4 807 m³ mr of precious species (15.18 % of national production), with a value of about $ 84 091 064 and $ 17 306 518, respectively, which places it with a contribution to GDP below 0.2 % registered at national level, according to the Ministry of Environment and Natural Resources (^{Semarnat, 2017}). Bursera simaruba (L.) Sarg. (Chakáh) and Metopium brownei (Jacq.) Urb. (Chechem negro) are common tropical species of commercial interest for the forestry industry (^{Silva et al., 2011}; ^{Gallegos et al., 2012}), they are widely distributed in the Peninsula of Yucatán (^{Herbario CICY, 2010}), are used for agroforestry for their growth and regeneration (Bursera simaruba) as well as for the production of furniture and handicrafts (Metopium brownei) (^{Román et al., 2016}; ^{ITTO, 2020}).

The accurate estimation of the stem volume, of the total volume when including the branches, and of the commercial volume based on the dimensions that the market demands, is important to develop the silvicultural practices proposed in the forest management plan according to the data of the forest inventories carried out (^{Corral-Rivas and Návar-Cháidez, 2009}; ^{Diéguez-Aranda et al., 2009}; ^{Gómez-García et al., 2016}; ^{Tamarit et al., 2017}; ^{Niño et al., 2018}).

Mathematical modeling is a support tool for forest resource managers, for which the adjustment of models that estimate the stem, total and commercial volume; o model, with precision, the ratio between the commercial volume and the stem (Vc / Vf), as well as its dynamics of decrease in the accumulated volume of a tree with the increase in height and the reduction of the diameter on the stem until reaching to the total height (volume ratio functions), in addition to describing the taper of the trees (diameter at different heights) by species, contributes to making management plans according to the specific and current conditions of the resources (^{Burkhart, 1977}; ^{Tapia and Návar, 2011}; ^{Barrios et al., 2014}; ^{Özçelik and Göçeri, 2015}; ^{Hernández-Ramos et al., 2018a}).

In recent years in Mexico, modeling has been carried out through functions that consider diameters at different heights on the stem and compatible commercial taper-volume systems for different species of the Pinus genus (^{Tapia and Návar, 2011}; ^{Quiñonez-Barraza et al., 2014}; ^{Tamarit et al., 2016}; ^{Silva et al., 2018}), or by volume ratio models for species such as Eucalyptus urophylla Blake in Huimanguillo, Tabasco, Swietenia macrophylla King in Quintana Roo and Pinus pseudostrobus Lindl. in San Juan Nuevo Michoacán, among others (^{Hernández-Ramos et al., 2018a}; ^{García-Espinoza et al. 2018}).

Commercial volume systems composed of ratio models are a tool to accurately calculate the commercial volume and the distribution of timber products from standing trees (^{Prodan et al., 1997}; ^{Hernández-Ramos et al., 2018a}; ^{García-Espinoza et al., 2018}). The estimation of these volumes is fundamental in the projection of the use of the timber resource and the yield for the forest industry, because they make it possible to estimate the financial income of the forest from forest inventory data, mainly the distribution of products per unit area. (^{Niño et al., 2018}; ^{Silva et al., 2018}); and they also contribute to avoiding the overuse of species of commercial importance in general, not only of tropical forests.

For the area of interest, ^{Hernández-Ramos et al. (2018a} and ^2018b) adjusted total volume models and generated a commercial volume system from ratio models, in addition to a segmented taper system (d) compatible with commercial volume (Vc) for Swietenia macrophylla King (mahogany) in tropical forests of Quintana Roo. ^{López et al. (2019)} adjusted a compatible system of commercial volume equations consisting of stem, branches and stem and thick branches of trees for 11 important tropical species in southeastern Mexico.

Traditionally, the equations to design these commercial volume systems have been adjusted using linear and non-linear regression models that assume normality, equality of variances and independence of errors (^{Carrero et al., 2008}).

However, obtaining reliable equations for some variables is complex, due to the variability that exists in these relationships due to the effect of the different forms of the stems and the conditions in which each individual develops (^{Quiñonez et al., 2012}). At present, the adjustment under the approach of mixed effects models (MEM) has become a viable option for this type of study (^{Corral-Rivas et al. 2014}), since with its incorporation in its parameters it allows grouping the sources of variability and reduce the variance of the error (^{Seoane, 2014}). This MEM technique includes in its formulation fixed parameters common to the entire population and specific random parameters for each level of grouping (^{Seoane, 2014}). As a result, more efficient, precise and reliable estimators result than the fixed parameters of the model and, therefore, estimates can be made with greater certainty of the variable of interest (^{Carrero et al., 2008}). For all the above, the objective of the present study was to fit a non-linear function with random effects that estimates the commercial volume for Bursera simaruba and Metopium brownei in central and southern Quintana Roo, Mexico.

Materials and Methods

The study was carried out in different forest areas of Quintana Roo, Mexico, which are located in areas with topoforms of low hills and rocky plains at an altitude of 20-180 m. The predominant climate is warm subhumid (Aw), with an average annual temperature of 22 °C and an average rainfall of 1 300 mm (^{Semarnat, 2014}).

From tours in the forest areas of the state, trees larger than 7.5 cm in normal diameter were selected randomly located in undisturbed areas of the tropical forest, where 188 specimens of Bursera simaruba and 133 of Metopium brownei were sampled, with which they were covered all the diameter categories found in these forest ecosystems (10-50 and 10-60 cm, respectively). Information on stump diameter (dt), normal diameter (d), total height (H), and partial diameters and heights on the stem (d _i and h _i ) were taken from each of them. The sample was chosen with monopod trees, without the presence of mechanical damage or pests or diseases, in addition to being at least 20 m away from roads, clearings or places with clearings to avoid edge effects. The descriptive statistics of the analyzed variables are arranged in Table 1.

The measurements of the diameters (di) at different heights on the stem (h _i ) of 0.3 m, 0.6 m, 0.9 m, 1.3 m and 2.5 m, were made directly with a caliper graduated to the millimeter, while the rest of h _i and d _i with a Bitterlich Telerelascope along the shaft where possible from 2.5 m to full height (H) or zero diameter (d ₀ ). For the calculation of the stem volume (Vf) the volumes of the different sections (logs) were added, which were determined by the Smalian formula and for the volume of the tip or last log (Vp) with the cone formula; the volume ratio (R) resulted from the ratio between the volume at different sections on the stem (Vi) and the Vf, for which the expression Vi/ Vf was used (^{Balzarini et al., 2005}; ^{Hernández-Ramos et al. , 2018a}; ^{García-Espinoza et al. 2018}).

As indicated in the Official Mexican Standard NOM-152-SEMARNAT-2006, a non-linear function (general model) for each species was adjusted, which is made up of two components: the first one includes a mathematical structure that corresponds to an equation of Vf (a ₀ +a ₁ (d ² H/1000) and the second one is an equation referring to the ratio between the volumes of the sections and the R (e ^{(-(a2)*((di/d)/1000)} e ^(a3*(di/d)) ), while these equations, when multiplied, estimate the commercial volume (Vc) (^{Schabenberger and Pierce, 2002}; ^{Balzarini et al., 2005}) (Table 2).The latter is defined by the dimensions that the sawmill industry or timber products demand, in which there is an interest in estimating, whether long dimensions or primary, short or secondary and waste or firewood. For the two species, the general model and the random effect additive to each combination of the fixed parameters were tested and eight expressions were obtained from the general model.

The general model (Table 2) was first adjusted by ordinary least squares (OLS) and two variants of ^{Spurr's (1952)} shaft volume equation: a linear one of the y =a ₀ +a ₁ X form, and the second one of the exponential type y = aX ^b form. In them, y is the response variable (Vf = stem volume), X represents the independent variable included in the model (d ² H), a and b correspond to the parameters to be estimated.

The incorporation of the random effect to the fixed parameters (a _i + u) was done per tree to explain the variability in its tapers, considering the biological meaning of each one; for example, by including the random effect in the parameter referring to the shape of the tree, the rate of change of the geometric body in the stem, the ratio of proportion between the stem volume and the commercial volume or in all three simultaneously within each of the components in the proposed non-linear function, where the probability value of the effect must be significant (P ≥ 0.05) (Table 1, 1-16).

The adjustment was carried out in the SAS 9.3^® program using the nlmixed procedure and the FIRO first order method, which defines with greater precision the approximation by likelihood of the values in the estimated parameters (^{Beal and Sheiner, 1982}; ^{SAS, 2011}) The distribution for both the dependent variables and the random effects (a _i + u) was specified in the programming as NORMAL to be adjusted by the Dual Quasi-Newton optimization technique and the Adaptive Gaussian Square integration method of the Statistical Analysis System (^{SAS, 2011}).

The selection of the best combination when including the random effects within the non-linear model, was made through the significance of the estimated parameters in both fixed and random effects at a value of 95 % reliability (P ≥ 0.05), and the best values of Verisimilitude-2Log and Information Criterion of Corrected Akaike (CIAC, for its acronym in Spanish) (^{García and Rapelli, 2011}; ^{Gómez-García et al., 2016}). Furthermore, for the best expression, the trend in the estimates was graphically verified by using the model structures separately in the stem volume (Vf = first mathematical structure) and volume ratio (R = second mathematical structure), and so joint to estimate the commercial volume (Vc) of each species or up to a diameter of 20 cm (Table 1, general model).

In the selected model, the homoscedasticity of the residuals was verified using White's test, under the hypothesis that the residuals are homogeneously distributed around zero (null hypothesis: Ho), while the alternative hypothesis refers to a homoscedastic distribution of the residuals (alternative hypothesis: Ha) at a 95 % reliability level (^{Martínez et al., 2006}). The quality of the estimates for the variables Vf, R and Vc was determined through bias (E) (^{Quevedo et al., 2010}) and the coefficient of determination (R ² ) and through the graphic analysis of the observed values against the values estimates of each variable (^{Pece, 1994}).

Results

When adjusting the general model and the variations in the inclusion of random effects, it was observed that both the general structure of the model and expressions 1, 2, 4 and 5 for Bursera simaruba and 1, 2, 4 and 5 for Metopium brownei, in some of its fixed parameters, they were not significant at 95 % reliability (P ≥ 0.05). While, in the combinations of inclusion of the random effects of expressions 6 and 8 in both species, the parameters related to mixed affects could not be estimated even when they tried to adjust by the Gauss-Hermite method (^{Pinheiro and Bates, 1995}). Therefore, these models were discarded for further analysis and were not included in Table 3. Due to this situation, only the combinations of inclusion of the random parameters of expressions 3 and 7 for Bursera simaruba and Metopium brownei trees are presented for subsequent analyzes, in which all were significant (P ≥ 0.05) when adjusted with 1 476 and 1 130 pairs of data, respectively.

For the two species, the combination in the mathematical structure that contains the Spurr stem volume model in its exponential form, presents the highest statistical values in the likelihood indicators -2Log and CIAC, for which the structures of model 7 were selected for Bursera simaruba and Metopium brownei, respectively. When contrasting these values with those obtained when adjusting the general model by ordinary least squares (OLS), an average statistical gain is observed when applying the mixed effects models approach of 5.6 % and 6.9 % for Bursera simaruba and Metopium brownei, respectively (Table 4). In addition, White's tests for both species were not significant (a < 0.05), which indicates that the residuals of both models are distributed in a homoscedastic way, which is why the null hypothesis (Ho) is accepted and Ha is rejected.

The biases estimated for the two species show good precision in the estimates when evaluating them individually for Vf and R, and as a complete model to estimate Vc, since in all cases they are less than 0.03 m³ (Table 5). In addition, when performing the graphical evaluation and by means of a regression of the observed values against those predicted by the variable of interest, no obvious problems were observed in the estimates, a situation that is reflected in the value of the determination coefficient (R ²) that in all cases was greater than 0.9 (Figure 1).

Figure 1 Trend of stem volume estimates, volume ratio and commercial volume for Bursera simaruba (L.) Sarg. (1a, 1b and 1c) and Metopium brownei (Jacq.) Urb. (1d, 1e and 1f) species. 

To exemplify the use of the structure that best fit the observed data of the commercial volume model, the values of normal diameter (d) = 30 cm and total height (H) = 15 m for Bursera simaruba, and d = 40 cm and H = 20 m for Metopium brownei, using model 7 with its parameters for each species, respectively, to estimate the stem volume (Vf) and volume ratio (R) when using the two sections of the model separately and when applying the complete model adjusted for commercial volume (Vc). From the general expression of the model, the first component corresponding to Vf is taken, which, when substituting for Bursera simaruba there is:

Vf = 0.0571 ((30 ² · 15)/1000)^0.9440 = 0.6663 m³

And for Metopium brownei: Vf = 0.0472 ((40 ² · 20)/1000)^0.9524 = 1.5049 m³.

Next, the minimum diameter for primary products is established, which in this case was 20 cm; then the second section of the general model is taken and the values are substituted, for which the following results for Bursera simaruba:

R =Exp(-9.1855·((20/30)/1000)·Exp(5.321*(20/30)))= 0.8085, y para Metopium brownei: R = Exp(-3.0967·((20/40)/1000)·Exp(6.3067*(20/30))) = 0.9644.

Finally, to obtain the Vc at the established minimum cutting diameter, the values of Vf and R are multiplied, or the proposed model of each species is applied directly,: for Bursera simaruba as follows:

Vc = [0.0571 ((30 ² · 15)/1000)^0.9440]([Exp(-9.1855·((20/30)/1000)·Exp(5.321*(20/30)))] = 0.5387 m³, y para M. brownei:

Vc = [0.0472((40 ² · 20)/1000)^0.9524] [Exp(-3.0967·((20/40)/1000)·Exp(6.3067*(20/30)))] = 1.4514 m³ .

Discussion

The MEM allows correcting the structure of variances-covariances associated with the classification variables of the trees (^{Littell et al., 2006}), since by grouping the information in class, under the MEM approach, the error is reduced, since the variability between the classes is explained by the inclusion of random parameters, which are estimated simultaneously with the fixed parameters (^{Castedo et al., 2006}; ^{Corral et al. 2019}). In this analysis technique (MEM), it assumes that the variance is specific for each classification level (^{Seoane, 2014}), so it is possible to increase the R ² and improve other goodness-of-fit criteria, such as the AIC and BIC.

In the evaluated models, random effects were included to account for the effect of classification covariates during adjustment in the equations. In this case, the random effect was specific at each level of the classification factor, which according to ^{De los Santos et al. (2006)}, ^{Tamarit-Urias et al. (2014)} and ^{García-Espinoza et al. (2019)} allows obtaining more efficient, precise and reliable estimators of the fixed parameters of the model and predicting specific random parameters of each sample unit, which reflect the deviation pattern with respect to the mean, but after adjustment, the values obtained from the fit for the random parameters are not used to construct the predictions, so the parsimony of the models is not modified.

In the case of the structure of the volume model 3 for Bursera simaruba, and 7 for M. brownei, all the parameters were significant (P ≥ 0.05) (Table 3), similar to that reported by ^{Quevedo et al. (2010)}, these authors analyzed the stem volume models under the mixed effects approach for Tabebuia rosea (Bertol.) In Venezuela, where a grouping by diameter category and transformation of the dependent variable was performed. In a similar way to the study carried out by ^{Gómez-García et al. (2016)}, who developed a profile function through mixed effects for Pinus sylvestris L. in Turkey, and where it was adjusted by the FIRO (first-order) method, of the R-project program, which is the same method used in this work.

For the choice of model 3 and 7 of Table 1 for Bursera simaruba and Metopium brownei, respectively, the description by ^{García and Rapelli (2011)} was taken as a reference, who reported that the models with the values closest to zero in the indicators of Verisimilitude -2Log and CIAC have a higher likelihood in the fit. In addition, these values were consistent with the evaluation criteria made by ^{Gómez-García et al. (2016)}; and by ^{Quiñonez-Barraza et al. (2014)}, who adjusted compatible commercial taper-volume systems for the main Pinus species from Durango, Mexico.

The bias values in the models proposed by species and the determination coefficients (R ²) obtained by contrasting the observed data against those predicted by the variable of interest (Vf, R and Vc), coincide with that reported by ^{Pece (1994)}, who used a volume ratio model to estimate the commercial volume in Eucalyptus pellita F. Muell. in Brazil and obtain biases when calculated by diameter category between 0.06 and 0.11. Furthermore, it is similar to that determined by ^{Hernández-Ramos et al. (2017)}, by proposing volume ratio systems for Eucalyptus urophylla in Tabasco and obtaining accurate estimates with lower biases in the variables of Vt, R and Vc than those obtained in this work (< 0.0077) and R² when contrasting the estimated values versus observed greater than 0.94.

Likewise, the precision of the estimates obtained by this work in commercial volume (Vc) are similar to that recorded by ^{Hernández-Ramos et al. (2018b)} when using volume ratio systems in Swietenia macrophylla trees in Quintana Roo and obtaining biases for volume ratio structures based on normal diameter and total height of -0.0121 and 0.0053, respectively. However, they differ with ^{García-Espinoza et al. (2018)} by proposing structures based on the volume ratio for P. pseudostrobus in Michoacán, and registering good accuracy in the estimates and a lower bias in the Vc compared to those obtained in this work (< 0.001 m³).

The graphical evaluation of the trends in the estimates (Figure 1) agrees with that obtained by ^{Tapia and Návar (2011)} who adjusted models of stem volume and taper for P. pseudostrobus in Nuevo León;^{Tamarit-Urias et al. (2016)}, by modeling the stem profile of nine species of the Pinus genus, and obtaining projections according to the trends of the observed data. It also coincides with those reported by ^{Tamarit et al. (2017)} in the construction of a cubing system for individual trees for the Quercus genus in Puebla and those reported by ^{Silva et al. (2018)}, who constructed compatible taper and commercial volume functions for three Pinus species in Chihuahua.

These tools generated at work can be complemented with the equations proposed by ^{García et al. (2017)}, which propose to estimate the stem volume through the dimensions of the stump diameter in Bursera simaruba and Metopium brownei in the tropical forests of Quintana Roo, and can be taken as an update of the equations proposed by ^{De los Santos (1976)} and ^{Patiño et al. (1994)} for these two species.

Conclusions

The proposed volume models based on diameter at different heights above the stem are reliable to estimate the stem volume when using the first part of the equation, the volume ratio when using the second part, and the commercial volume when used together in these two species. This is a very significant advantage if it is compared with equations generated more than 30 years ago to estimate clean stem volume or total stem volume for large trees, which does not correspond to the diameters that are used today. With the mathematical structures proposed by species, it will be possible to build product distribution tables to project the timber yields of areas under forest use, according to the current conditions of the Bursera simaruba (Chakáh) and Metopium brownei (Black Chechem) species in tropical forests of the central and southern state of Quintana Roo, Mexico.

It is necessary to carry out this type of study in Quintana Roo, in order to cover at least the 30 main species that are most important from a commercial point of view.