Study area and dasometric data

The study area pertains to the municipality of Zacualpan, Veracruz, Mexico (20° 25’ and 20° 31’ N, 98° 24’ and 98° 25’ W, with 1670 m altitude). The climate is temperate humid with an average temperature of 18 °C and mean annual rainfall of 1900 mm; soil is Feozem.

In 2011 50 permanent observation plots were established, which were re-measured in 2012. Of these, in 2015 14 cleared plots were found, 26 with partial mortality (intra-specific competition), seven with total use and three with total death of the plants by dominance of weeds. The plots are square, 20x20 m (400 m²). In each period normal diameter of all of the trees was measured (DN in cm) and the height of the 5 dominant trees (A in m) and the age of the plantation was measured (E in years).

The state-variables by hectare were: average dominant height (A, m), total volume (V, m³), squared diameter (Dq, cm), basal area (AB, m²) and number of residual trees (NA), all referred to the site index (SI) as criteria for classifying the site quality (^{Clutter et al., 1983}).

Explicit growth system

The explicit SCRM was structured by integrating models of prediction and projection for the variables A, AB, V and NA). The results were models of estimating future yield per hectare and by station quality determined through a site index (^{Clutter et al., 1983}; ^{Burkhart and Tomé, 2012}). The construction of the system initiated with the selection of a generic model for dominant height and afterwards using the algebraic difference method (ADA) a future condition was derived as a function of an initial situation (^{Pienaar and Rheney, 1993}; ^{Jordan et al., 2006}). Dominant height was fitted with non-linear regression independently and the AB, V and NA as a compatible system (^{Santiago-García et al., 2013}).

Dominant height

In the fitting of the models in dominant height (A) and site index (SI) 76 pairs of non-overlapping data were used that correspond to the average of the dominant heights with their ages, which varied from 3 to 21 years. From the prediction model, the models of projection and site index were derived. The algebraic difference method (ADA) generates families of anamorphic and polymorphic curves of the form A ₂ = f (A ₁, E ₁, E ₂, β) where: A ₂ is the projected dominant height or site index at age E ₂ (projection age), A ₁ is dominant height (at a base age of 15 years) measured at the initial age (E ₁) and β is the vector of regression parameters (^{De los Santos-Posadas et al., 2006}; ^{Quiñonez-Barraza et al., 2015}). The expression of Schumacher-Korf was used in its form of guide curve (1) and polymorphic-2 (2) with growth hypothesis in the parameter α₂ of the generic model (1) (^{Schumacher, 1939}).

(1)

(2)

where A: dominant height of prediction, A₁: dominant height at the initial age, A ₂: dominant height at the projected age, α₀, α₁, α₂: regression parameters, E₁: Initial age or base age (15 years), E ₂: Projection age, e: exponential function and ln: Napierian logarithm.

Basal area

Twenty four permanent observation plots were segregated from the sample of 2015 due to thinning, definitive cuts or total death. A total of 73 pairs of non overlapping data were used for the fitting of the basal area, volume and number of trees per hectare. The proposed model of basal area (3) is a function dependent on the dominant height, age and number of trees per hectare. With this expression the model (4) was generated for estimating the basal area at the projected age (AB ₂) (^{Gonzalez-Benecke et al., 2012}). The expressions are as follows:

(3)

(4)

where AB₁, A₁, E₁, NA₂: basal area, dominant height, age and number of trees per ha at the initial age, AB₂, A₂, E₂, NA ₂: basal area, dominant height, age and number of trees per ha at the projected age and β₀, β₁, β₂, β₃: parameters to be estimated.

Volume

For the volume the pair of equations (5 and 6) was used which describe the volume of prediction and of projection, respectively:

(5)

(6)

where V₁, AB₁, A₁, are volume, basal area and dominant height in their initial condition, V₂, AB₂, A₂: volume, basal area and dominant height in the projected condition and γ₀, γ₁, γ₂: parameters to be estimated.

Number of residual trees

Outstanding in the model for the number of residual trees, is a parameter that represents the mortality rate (ω₁); this parameter multiplies at the age of reference. The survival can be expressed in terms of prediction (7) and projection (8) in the following form:

(7)

(8)

where NA₁ , E₁, NA₂, E ₂: Number of trees per ha and age in the initial and projected condition, respectively and ω₁, ω₂: parameters of the regression.

System of implicit growth using the fdp Weibull

^{Weibull (1951)} published a probability density function (fdp) and its power is such that it has even been widely applied in dasonomy (^{Villaseñor-Alba and Villanueva-Morales, 2000}). To fit the Weibull function and model the distribution of diameter classes there are methods such as that of linear and non-linear regression, percentiles and maximum likelihood (^{Torres-Rojo et al., 1992}; ^{Torres-Rojo et al., 2000}). The method of estimation of parameters through percentiles has shown good results (^{Santiago-García et al., 2014}).

The state-variables proposed were: basal area per hectare (AB, m² ha^-1) and mean squared diameter (Dq, cm) given by:

(9)

The height of the unmeasured trees was obtained indirectly with the equation H _i =0.880151xA ^1.026997x(1-e ^{(-1.075746xDN)} ) that relates dominant height (A) and normal diameter (DN) to estimate the individual total height. Total tree height with bark was obtained with the Schumacher model: Vi=0.000114×DN ^1.718158 ×H ^0.939592 (^{Uranga-Valencia et al., 2015}) where H=total tree height and DN its normal diameter. The proposed percentiles were as follows:

(10)

(11)

(12)

where p0, p63, p95 correspond to the zero percentile, 63 and 95, D_q: mean squared diameter, AB: basal area and A: dominant height, φ₀, φ₁, φ₂, ω₀, ω₁, ξ₀, ξ₁, ξ₂ are the parameters to be estimated in the percentiles system. The Weibull 3P probabilities distribution function is given by the probabilistic function:

(13)

(a≤x≤∞); =0, otherwise, where a is the parameter of location, b of scale and c of form. The function defines the value of probability density associated to a random variable, which for this study was normal diameter (x=DN) (^{Santiago-García et al., 2014}). The accumulated distribution in its closed form is:

(14)

(a≤x≤∞); =0, otherwise.

The values of the parameters b and c are positive, whereas a can be negative, zero or positive. However, in the modeling of distribution of diameter classes, a should not be negative (^{Clutter et al., 1983}). Knowing the parameters, the density function for a given interval can be obtained with:

(15)

where P is the proportion of trees of the assigned diametric class, L and U are the lower and upper limit of the diameter class respectively, e is the exponential function and X is the center of the diametric category.

The parameters were estimated by the moments method (^{Pienaar and Rheney, 1993}), which has been tested with satisfactory results (^{Santiago-García et al., 2014}; ^{Quiñonez-Barraza et al., 2015b}), utilizing the following equations:

Parameter a of location was determined with:

(16)

Parameter c is given by c=ln [((-ln(1-0.95))/(-ln(1-0.63)))/ln(p95-a)/p63-a))]. For the parameter b the following expression of the second moment of the Weibull distribution was used:

(1 (17)

where Γ₁= Γ(1+1/c), Γ₂= Γ(1+2/c), Γ(.) is the Gamma function, the other variables were previously defined already.

Variables of state by diametric category

The number of trees per diametric class was obtained as follows: NA _(CD) =P(L<X<U)xNAHA is the function of density of probability for the interval of this same category (CD) and NAHA: number of trees per hectare residuals at a determined age. The volume per hectare per diametric category was obtained with the expression V _(CD) =NA _(CD) xVi, where NA _(CD) is the number of trees of the distribution of the diameter class.

Technique and test of goodness of fit

The fitting of the dominant height model and site index was done using non-linear progression with the MODEL procedure of the software SAS/ETS^® (^{SAS, 2014}). The coefficient of determination fitted by the number of parameters (R² _adj), the sum of squares of the error (SCE), the mean square of the error (CME), the level of significance of the parameters (Pr>|t|), and the logical modeling of the curves with respect to the graph of data observed in the field were the criteria of goodness to validate the model (^{De los Santos-Posadas et al., 2006}).

The compatible models of basal area, volume and mortality were fitted simultaneously with seemingly unrelated regression (SUR); this technique assumes that the variables are independent from each other. The test of Kolmogorov-Smirnov (KS) was used to qualify the goodness of fit of the percentiles. To discard plots a significance level of α=0.05 and 0.1 was used (^{Santiago-García et al., 2014}; ^{Quiñonez et al., 2015a}).

Dominant height and explicit yield

For dominant height the models of Hossfeld IV, Schumacher-Korf and Chapman-Richards have the best fits and from them are derived families of site index (SI) (^{De los Santos-Posadas et al., 2006}; ^{Quiñonez-Barraza et al., 2015}). Model (2) fitted to predict dominant height coherently describes the growth pattern of the permanent observation plots (Figure 1). This made it possible to derive four site qualities labeled by the site index (SI, m) at a base age of 15 years, with SI=10 for low productivity stands and SI=25 for stands with higher production potential. SI=20 was considered the most representative condition. The estimator of the parameter α₀ related to the asymptote has a rejection value for the model (1); this would imply that the height of the trees would reach 125 m. In contrast, model (2) has a projected value of 30.8 m height (Table 1).

Figure 1 Family of growth curves for dominant height by site quality with polymorphic pattern vs the growth pattern of the data observed in permanent plots. 

Table 1 Regression parameters and tests of goodness of fit of the curves of Schumacher-Korf (1) and of the polymorphic curve-2 for dominant height. 

SCE: sum of squares of the error. CME: mean square of the error; α₀, α₁ α₂: parameters estimated by non-linear regression.

Basal area, volume and number of residual trees

^{Pienaar and Rheney (1993)}, ^{Galán et al. (2008)} and ^{Magaña et al. (2008)} generated systems of knowledge and timber yield (SCRM) relating models of dominant height (AB, m² ha-¹), volume (V, m³ ha^-1) and mortality or number of trees (NA) with different statistical techniques, with good results. The system of equations of the present study generated adequate results. The criteria of goodness of fit (Table 2) are acceptable and with high levels of significance in their parameters (Table 3). With the system of equations and the estimated parameters in combination with the family of polymorphic curves, the table of timber yield was constructed for station quality and average quality (Table 4). The value of R² of the volume model is low, which implies a very high variation among the values of initial density. At year-one the parameter α₀ was fixed at 1100 trees per hectare, which is the initial density of most of the plantations of P. patula in Zacualpan, Veracruz. The projection of the growth curves congruently describes the pattern observed in the permanent observation plots (Figure 2A).

Table 2 Goodness of fit statistics under seemingly unrelated regression (SUR) of the system of equations of prediction-projection for basal area, volume and mortality. 

Table 3 Estimated parameters and significance level of the parameters with seemingly unrelated regression (SUR) of the prediction-projection system of equations for basal area, volume and mortality. 

Pr>|t|: test of significance of parameters; <0.001 is highly significant. β₀, β₁, β₂, γ₀, γ₁, ω₀, ω₁. Parameters of regression.

Table 4 Table of growth and yield AB, V, and NA for site quality, generated by the implicit method and non-linear seemingly unrelated regression. 

NA ^†: number of trees per hectare, SI: site index in meters at the base age of 15 years and age in years.

Figure 2 (A) Explicit modeling of growth in volume for different site quality superimposed over the trajectory of the data of re-measured plots. (B) Mean Annual Increment and Current Annual Increment (IMA and ICA) for site quality. 

The technical turn for site quality with SI 10, 15, 20 and 25 occurs at 21, 18, 14 and 11 years, respectively. The mean annual increment for these qualities is of 7, 12, 18 and 26 m³ ha^-1 year^-1 (Figure 2B).

Having three measurements of the same sampling unit and having plots with representative ages in the plantations resulted in a good fit of the models. Consequently it made it possible to adequately describe the growth pattern of the variables of interest.

In forests of the municipality of Zacualtipán, Hidalgo, regenerated naturally with P. patula, at the age of 20 years a yield in volume of 10.2 m³ ha^-1 year^-1 has been registered in site qualities of SI 26 m at a base age of 40 years (^{Santiago-García et al., 2013}). The above means that the plantations show a better yield when compared with natural forests.

The rate (parameter) that translates in mortality of the number of trees per hectare is the result of the intra-specific competition of the stand. In Zacualtipán a parameter of -0.0337 (^{Santiago-García et al., 2013}) was reported and in our investigation it was a little higher (-0.04). Growth in basal area in combination with a thinning sequel is useful for simulating models of response to thinning.

System of growth and implicit yield with the fdp Weibull

The statistics of goodness of fit and the estimated values in each parameter make it possible to apply the proposed models. The KS test with α=0.05 [C(α)=1.36] indicated the rejection of four of the 50 plots (8 %) and with α=0.1 [C(α)=1.22] a rejection of five (10 %). These results are similar to those of ^{Quiñonez et al. (2015a)} and ^{Santiago-García et al. (2013)}, the latter with 1.58 and 3.96 % of rejection for α=0.05 and α=0.1, respectively. The KS test confirms that the data of plots has a Weibull distribution. To this respect, at a higher level of significance the contrast of the deviations is reduced (^{Santiago-García et al., 2014}).

The SUR procedure utilizes the correlations of the errors of the models to improve the fit. To model Weibull 3P at least two percentiles are required, and their selection depends on the sample (^{Borders et al., 1987}; ^{Santiago-García et al., 2014}). The efficiency of the percentile 63 and 95 is shown in Tables 5 and 6. By being centered the bias of prediction is reduced. Similar studies have used percentiles 50 and 90, 40 and 24 and 93 (^{Santiago-García et al., 2014}).

Table 5 Statistics of goodness of fit of the models of percentiles for the Weibull function. 

p0, p63, p95: percentile zero, 63 y 95.

Table 6 Statistical parameters of the models of percentiles for the Weibull function. 

φ₀, φ₁, φ₂, ω-1 , ω-2 , ξ₀, ξ₁, ξ₂: Parameters estimated through non-linear regression.

The diameter distribution at the age of 20 years indicates that stands with average site quality (SI 20) will have 23 % more trees for saw-timber (DN>25 cm) compared to sites of poor quality (SI 10), which is reflected in the value of production, the larger the diameter, the higher the sale price of the product (Figures 3 and 4). In SI 25 at the age of 10 years it is projected to obtain 51 trees with DN<30 cm; for poor site class zero trees are projected with this same diameter class.

Figure 3 Diameter distributions estimated at different ages of the stand by station quality with the Weibull probability distribution function. SI: 10, 15, 20 and 25 m, Age 5, 10, 15, 20, 25 and 30: is the plantation age in years, NA ha^-1: number of live trees per hectare, diametric category in intervals of 5 cm. 

Figure 4 Timber yield (m³ ha^-1) by site quality: SI 10, 15, 20 and 25, against the age of the plantation. <15: total volume of wood in trees less than 15 cm of DN. >30 : total volume in trees with DN higher than 30 cm. Sum: accumulated volume integrating all of the diameter classes  

Promising results are expected for SI of 20 m; because at the age of 15 years 429 trees ha^-1 are predicted with DN<30 cm and 151>30 cm ha^-1, which corresponds to a volume of 155 and 103 m³ ha^-1, respectively. ^{Santiago-García et al. (2014)} present for the same species (average site) but with natural forest management of Zacualtipán, Hidalgo, a projection at 20 years of 70 trees with DN<30 cm and 750 individuals with DN<30 cm.

Results show a wide variation in the state-variables, which is motivated by site quality. In SI 10 the technical rotation is projected at 21 years with similar timber production.