Introduction

The classification of the productivity of forest areas is the basis for planning the establishment and management of commercial forest plantations (PFC, for its acronym in Spanish) (^{Carrero et al., 2008}), since the difference in productivity within the localities determines critical aspects of management, such as age (E) to apply forestry treatments, harvest shift, harvesting techniques or profitability indicators (^{Andenmatten and Letourneau, 1998}; ^{Gaillard et al., 2002}). Given this scenario, it is essential to know the growth rate of the species planted for market purposes (^{Aguirre-Calderón, 2015}).

The most widely used methodology to classify the station quality of forest stands is the site index (IS, for its acronym in Spanish) (^{Álvarez et al., 2004}; ^{Cornejo et al., 2005}; ^{Madrigal et al., 2005}; ^{Hernández-Ramos et al., 2014}), defined as the dominant height (Ad) that the stand reaches at a base (E ₀ ) or reference age. Ad in a monospecific and contemporary stand is little affected by density, follows a certain pattern and has a high correlation with the production of timber volume (^{Clutter et al., 1983}; ^{Torres and Magaña, 2001}).

To represent the IS, Ad curves are used that express the average growth pattern of the trees through non-linear equations, from which season quality or productivity curves are constructed, which can be proportional to each other (anamorphic) or possess particular growth patterns for each IS (polymorphic) (^{Clutter et al., 1983}; ^{Torres and Magaña, 2001}; ^{Hernández-Ramos et al., 2015}).

Until 2014, there was an area of 2 722 ha in the state of Quintana Roo covered with PFC, where Cedrela odorata L. was established in 2 155 ha as a priority species (^{Conafor, 2016}). The lack of information on the productivity of the plantations and on biometric tools for their management, has caused that the maximum potential yield and profitability are not reached. Therefore, and due to the importance of the calculation of the IS for the planning of forest management of the wooded areas under cultivation, and from the hypothesis that there is a differential in the growth potential between localities and given that equations can be adjusted that model the functional relationships between the Ad and the E of the trees, the objective of adjusting and comparing IS equations for PFC of C. odorata in Quintana Roo, Mexico was stated.

Materials and Methods

The study was carried out with data from the C. odorata PFCs located in the municipalities of Bacalar, Felipe Carrillo Puerto, José María Morelos and Othón P. Blanco, Quintana Roo, Mexico. The predominant climate is warm subhumid with rains in summer and part of winter, the thermal oscillation is less than 5 °C (Aw (x') i (^{García, 1983}). The average annual temperature is 25.5 °C and the average precipitation annual of 1 260 mm, with a pronounced dry season between December and May (^{INEGI, 2016}).

The PFCs are located in semi-flat lands, in soils that according to the Mayan terminology and its equivalent according to the IUSS Working Group ^{WRB (2007)} correspond to: Tzek’el (Leptosol), Kakab (Rendzina); K’ankab (Chromic Luvisol), Yax-hoom (Pelic Vertisol) and Ak’alche (Vertisol gleyco).

24 PFCs from one to 28 years of age were evaluated (Figure 1) using 400 m² temporary plots, in which the dominant trees in height were selected and 1 941 pairs of Ad-E data were obtained. Individuals without damage caused by pests, diseases or mechanical effects were chosen. A Sunnto^TM clinometer was used to measure Ad, and age was calculated based on the owners' information on the year of establishment.

Figure 1 Location of the commercial forest plantations of Cedrela odorata L. in the study area. 

The growth models of Schumacher (1), Champan-Richards (2) and Weibull (3) reported in the specialized literature (^{Kiviste et al., 2002}; ^{Burkhart and Tomé, 2012}; ^{Panik, 2013}) were adjusted to the refined database and which have shown satisfactory results in some studies to describe growth in Ad-E and IS (^{García et al., 1992}; ^{García et al., 2007}; ^{Hernández-Ramos et al., 2014}; ^{Murillo-Brito et al., 2017})(Table 1).

The information was processed in the SAS 9.4^® statistical package with the Model procedure and the maximum likelihood technique (^{SAS, 2016}). The determination of the goodness of fit of the models was based on the results of quantitative analysis, used in modeling studies and forest biometrics (^{Corral-Riva et al., 2007}). The lowest value of the root of the mean square of the error (REMC) and the highest value in the adjusted determination coefficient (R ² _aj. ) Were considered; while, the average bias indicators (E-) and the aggregate deviation in percentage (DA %) were included, in addition to the significance of the parameters and the distribution of residuals (^{Prodan et al., 1997}; ^{Lencinas and Mohr-Bell, 2007}; ^{Hernández-Ramos et al., 2020}). To verify the assumptions of normality of the data and homoscedasticity of the residuals, the Shapiro-Wilk (SW) test and the graphical distribution of residuals were used, respectively (^{Hernández-Ramos et al., 2020}).

Since the quality of fit does not necessarily reflect the adherence of the estimates to the observed data (^{Kozak and Kozak, 2003}; ^{Quevedo et al., 2004}; ^{Hernández-Ramos et al., 2015}), to reinforce the selection criteria of the best equation, a graphical comparison of the guide curves superimposed on the observed data was made and the one with the best fidelity was chosen. In addition, the model that showed a sigmoid trend was selected, a point of inflection in which the current (ICA) and average (IMA) increases become equal, since some do not usually adequately describe the growth of the trees and an asymptote or maximum possible growth at an infinite age. The crossing point at which ICA = IMA, was obtained according to the methodology of ^{Kiviste et al. (2002)}.

With the base equation selected for its adjustment and prediction, the guide curve in Ad was generated, considering that when E is equal to the base age (E₀), then Ad is equal to IS. From the IS equation, the parameter referring to the asymptote was cleared and the expression was replaced in the basic equation, with which the asymptotic parameter is considered implicit and the shape parameters were common for all sites (^{Clutter et al., 1983}); thus, by varying the E and the IS classes, and keeping the E ₀ constant, the anamorphic curves were generated. For the polymorphic curves, some of the parameters of the exponents were cleared from the IS equations, either the one associated with the growth rate or the site and they were substituted in the base equation. In a similar way, by varying the E and IS classes, and keeping E ₀ constant, the polymorphic curves were generated (^{García et al., 2007}).

By establishing an E ₀ equivalent to the shift established in 25 years and according to the observed data dispersion interval, IS classes were assigned for the station quality curves, while the values of the IS classes were substituted in the model and E ₀ (^{Clutter et al., 1983}; ^{Carrero et al., 2008}; ^{Tamarit-Urías et al., 2014}; ^{García-Espinoza et al., 2019}).

Results and Discussion

The dispersal trend of the 1 943 pairs of data used Ad-E of C. odorata trees show a sigmoid trend, where an increase in the growth of Ad in the PFCs can be observed as age increases (Figure 2).

Figure 2 Dispersion of the observed dominant height-age data of trees established in commercial forest plantations of Cedrela odorata L. in Quintana Roo, Mexico. 

The age interval ranged from one to 28 years and the Ad from 0.49 to 17.66 m, with minimum and maximum values between 0.22 and 19.35 m, respectively. The variation in Ad of the trees for each age is low, since they fluctuate around 0.71 % to 42.26 %, but most of the plantations vary between 2 % to 11 % (Table 2).

Table 3 shows the goodness of fit criteria and in which models (3) and (2) reached the highest values in the adjusted coefficient of determination (R ² _adj .). The lowest RCME and all their parameters were significant ( α <0.05). In this regard, ^{Gujarati (2010)} argues that values higher than 0.7 and 0.8 are acceptable, so it can be said that the three models explain in a good way the variability of the data.

The values in the determination coefficients are similar to those found by ^{García et al. (1998)} for Swietenia macrophylla King, with a R ² _aj value equal to 0.99 when using the same models, but higher than that calculated by ^{García et al. (2007)} of 0.84 in this criterion for C. odorata when using the Schumacher model as the one that best explains the variability of the information, as indicated by both authors when modeling the SI in PFC of the state of Quintana Roo, Mexico.

Several investigations record interesting results. Thus, ^{Minoche et al. (2017)} obtained different values when adjusting five models for dominant height and generating anamorphic IS curves using the algebraic difference technique in Tectona grandis L. plantations in the state of Tabasco, Mexico, where the best adjusted model was that of Sloboda, when determining an R ² of 0.998. ^{Galán et al. (2008)} adjusted and established an R ² of 0.870 and 919 in the Chapman-Richards model to predict and explain the dominant height growth and generate anamorphic IS curves of C. odorata and Tabebuia donnell-smithii Rose at a base age of 10 years that would cover the age ranges of 3 to 10 years. ^{Pérez-González et al. (2012)}, when adjusting seven IS models for S. macrophylla plantations of seven to 16 years in Tabasco, Mexico, concluded that the Korf model was statistically better to estimate growth in Ad at the E ₀ of 10 years in trees of the studied species (R ² = 0.80). ^{Martínez-Zurimendi et al. (2015)} calculated a value of 0.58 for the Chapman-Richards function and 0.59 for the Weibull model in Gmelina arborea Roxb. plantations. ^{Murillo-Brito et al. (2017)} reached a value of 0.50 in the Schumacher model for C. odorata plantations, both in Tabasco, Mexico.

Figure 3 shows the distribution of residuals, defined as the difference between the real observed value and the estimated value through the generated equation, in which the analysis shows that the three models have a good fit and a regular distribution around zero.

Figure 3 Distribution of residuals of growth models to generate site index for Cedrela odorata L. in Quintana Roo, Mexico. 

Model (3) had the highest adjustment because the deviation of the predicted and observed values (E-) is low, since the model overestimates by -0.0131 m at the individual tree level and at the sample level a DA =-0.2810 %. The negative values of the mean of the residuals are explained because there is a slight overestimation in the data of the predicted Ad as a function of the E for C. odorata plantations (Table 4).

When verifying the Shapiro-Wilk test to evaluate normality, the best models were (2) and (3), since they presented a value of SW = 0.89 close to 1. The Schumacher model obtained a lower value of SW = 0.88, but despite this, it can be said that it showed normality of the errors. These results are lower than those obtained by ^{García et al. (2007)}, obtaining a value of SW = 0.97 for C. odorata forest plantations in Quintana Roo, but despite having lower values, in the functions used the regression assumption is not being violated and therefore they are acceptable. The three models fulfilled the assumption of a homogeneous distribution of the residuals (Figure 3).

These results are higher than those of ^{Martínez-Zurimendi et al. (2015)} in the referred investigation with G. arborea in Tabasco, when they used the Chapman-Richards and Weibull models, in which the bias values were 4.05 m and 4.13 m, respectively; in the case of DA %, the numbers were higher too, 16.84 % for the Chapman - Richards model and 16.83 % for the Weibull model. In the aforementioned study carried out by ^{Murillo-Brito et al. (2017)} with the Schumacher model for forest plantations of C. odorata also in Tabasco, the figures were lower, the bias was -0.0398 m and the DA % of 1.9678 %, which resulted in a better fit than that obtained in this research with the same model. The equations obtained to estimate the guide curve were the following:

In equation (3), the parameter β ₀ represents the asymptotic value or maximum height that trees can reach (theoretically at an infinite age), in this case at 28 years of age, when they measure up to 16.63 m, since it represents the universe of Ad-E data used to fit the equation. For its part, the parameter β ₁ is considered a site-dependent parameter and represents the growth rate or the number of units that β₀ increases when we vary one unit in β ₁ , adjusting the scales by the exponential expression. The parameter β ₂ is associated with the site.

The guide curves and the observed values were plotted to analyze the fidelity of the models, which is observed in Figure 4a, where the guide curve of the model (3) exhibits a better trend since it tends to pass in the middle of the cloud of observations, so it more precisely describes the database used in its adjustment. In the case of model (1), the guide curve at younger ages tends to underestimate and at older ages to overestimate, while model (2) overestimates on a large scale at younger ages. In this sense, the research by ^{García et al. (2007)} for the same species indicated that the best model that showed an adequate trend in the guide curve was that of Schumacher, and in the study by ^{Murillo-Brito et al. (2017)}, the models that showed the best fidelity of the guide curves graphically with the observed data were Hosfeld IV, Monomolecular, Logistic and Gompertz.

Figure 4 Growth and increment of curves to generate the site index for Cedrela odorata L. in Quintana Roo, Mexico. 

According to the above statements, the best model is (3), since the guide curve represents the average of the observed data, it presents an inflection point where ICA = IMA, according to this it can be said that it adequately describes the growth in tree height as a function of age (Figure 4b). Model (1) also has a culmination point, but the trend of the guide curve is not suitable for this study. In the case of model (2), it does not present an inflection point, which means that it does not adequately describe the database used in its adjustment and cannot be used to carry out the IS curves.

In order to calculate the age at which ICA = AMI, the formula proposed by ^{Kiviste et al. (2002)} was used; therefore, a technical growth shift in Ad of eight years was obtained for this species (Figure 4b). This is similar to what was obtained by ^{García et al. (2007)} for C. odorata in the southern region of the state of Quintana Roo, defining a culmination of growth in Ad of 6.38 years. On the other hand, the results differ from those obtained by ^{García et al. (1998)} who determined five IS (12 m, 14 m, 16 m, 18 m and 20 m) at an E ₀ of 30 years in S. macrophylla plantations in Quintana Roo, with the Chapman and Richards and Weibull models, which in their polymorphic version overestimate height at young ages and underestimate it at older ages, however, this is possible due to the growth habits and ages of the plantations.

With equation (3), IS curves were constructed at an E ₀ of 25 years (Figure 5), three IS classes were used (14 m, 17 m and 20 m). The anamorphic (a) and polymorphic (b) expressions using the guide curve method are expressed below:

(a) Anamorphic: Ad=IS1-e-0.0352429(25)1.5524551-e-0.0352429(23)1.552455

(b) Polymorphic:      Ad=16.636711-e-23250.035244ln1-IS16.63671

Figure 5 Anamorphic (a) and polymorphic (b) growth curves of site index (IS, m) at the base age of 25 years in commercial forest plantations of Cedrela odorata L. in Quintana Roo, Mexico. 

In Figure 5, the anamorphic (5a) and polymorphic (5b) -type IS curves can be seen superimposed with the observed data. These curves make it possible to classify C. odorata plantations according to their productivity and thus clearly differentiate their seasonal quality.

The curves are different in the growth trend of the curves generated by ^{García et al. (1992)} for S. macrophylla, since they obtained five classes of IS curves with the Chapman-Richards model, and by ^{García et al. (2007)} for C. odorata in Quintana Roo with the Schumacher model, which has a different geometry. For their part, ^{Murillo-Brito et al. (2017)} for C. odorata in the state of Tabasco, generated three IS curves with the Hossfeld IV model, resulting in site indexes of 11 to 19 m.

Figure 5 indicates that the anamorphic curves cover the trend of the dispersion interval of the data, which shows that the proposed method is correct to satisfy the aim of this study, although the IS curve of 20 slightly overestimates at older ages. In the case of the polymorphic growth curves (Figure 5b), they were assigned other labels different from IS, since with the 14 m, 17 m and 20 m values of the curves, they tend to overestimate the observed data on a large scale. The values of the 14 m, 15 m and 16 m labels used to construct the polymorphic curves vastly underestimate the observed data from Ad-E. With these results it can be asserted that the curves that presented the best degree of association with the observed data were those of the anamorphic type.

From the guide curve, the E ₀ defined at 25 years and the three IS values, established at 14 m, 17 m and 20 m, it was possible to develop a system of anamorphic curves that represent the growth pattern for C. odorata in the central and southern region of Quintana Roo, Mexico.

With the expressions obtained, the IS of any plantation growing in the study area can be qualified. Assume a 23-year-old plantation with 14.55 m of average Ad and that the cutting shift in the region is 25 years, which will be defined as the reference age. These data are substituted in the expression and will indicate that the IS of that plantation is 16.57 m at the age of 25 years.

Conclusions

The Weibull model was the best fit for dominant height-age (Ad-E) data and for generating anamorphic and polymorphic site index curves of C. odorata in Quintana Roo, Mexico.

The polymorphic curves underestimate the growth in height of the dominant trees of the stands with the three site index labels (IS) used, while the anamorphic curves represent greater fidelity with the observed values, since they cover the trend of the interval of dispersion of data. For this reason, they are considered the most appropriate to qualify the IS in commercial forest plantations of C. odorata in Quintana Roo.

The information generated serves to form part of the planning methods for forest managers and technicians in forest management and exploitation activities, but the values obtained should be applied within the range of the variables Ad-E, which is from one to 28 years old.