1. INTRODUCTION
To know the tapering of tree shafts is important for planning and implementing forest activities, especially when one wants to classify the production through the wooden logs. By associating it with the growth factor, the mathematical description of the tapering allows to dynamically infer the quantity and dimension of the logs, simulating scenarios that contribute to studies on the economic viability of the forest (^{Costa et al., 2016}). Simulations of different scenarios are possible through application of linear and nonlinear modeling techniques. Thus, it is possible to quantify the multiple products from the settlements.
In this context, nonlinear models are preferable to describe biological phenomena, since they usually provide more precise estimates when compared to linear models. Besides, for some specific cases, this accuracy in the estimate is also associated to direct interpretation of the parameters, linked to the flexibility of application that they provide (^{Santos, 1996}). Some authors, such as ^{Mendonça et al. (2007}), ^{Pires & Calegario (2007}) and ^{Horle et al. (2010}) verified the superiority of nonlinear models regarding linear ones in the modeling of tree profiles.
A prediction model usually needs a large and representative data set, addressing different characteristics for the construction and adjustment stages. Therefore, traditional models with only fixed coefficients do not consider the possible variation of parameters between different groups. Thus, mixed models can deal with this variation, considering some model parameters to be random (^{Pinheiro & Bates, 2000}). Given this context, forest science already has major contributions in this topic, such as ^{Calegario (2002}) and ^{Calegario et al. (2005}) on de modeling of growth and forest production of Eucalyptus sp.; ^{Trincado et al. (2007}), ^{Meng et al. (2008}) and ^{Mendonça et al. (2015}) in the description of the hypsometric relation; ^{Vismara et al. (2016}) in volume prediction; and ^{Carvalho et al. (2014}), ^{Arias-Rodil et al. (2015}, ^{2016}) and ^{Môra (2015}) in the tapering modeling.
Therefore, with the advances in computing and biometric techniques, estimates can be obtained more accurately and more rapidly by generating the parameters of the models. Given this, the objective was to assess and compare the estimates of the partial volumes of Pinus taeda obtained by the nonlinear mixed modeling compared to the traditional one.
2. MATERIAL AND METHODS
2.1. Area of study
The study was performed in the municipality of Campo Belo do Sul-SC, in the region of the company Florestal Gateados Ltda. The place is circumscribed by the following geographical coordinates: S28°03’26” W50°46’13”. According to Koppen’s classification, the climate in the study area is predominantly Cfb, mesothermal, humid subtropical, with cool summers, without defined dry seasons, with occurrence of severe frosts and average annual temperature of 16 °C. Annual precipitation varies from 1,300 to 2,400 mm and the average altitude is 950 m. The most representative soil in the place is Haplic Nitisol, with associations of Cambisol and Litholic Neosol in the sloping areas (^{Embrapa, 1988}).
2.2. Data collection
The dendrometric variables for behavioral analysis of the trunk profile are distributed in a data set of 558 Pinus taeda trees in four forest sites provided by the company, ranging from 11 to 31 years-old.
The measured trees were selected based on the diametric distribution of the settlements from the forest inventory data. These selected individuals were felled and cut and their diameters were measured at different heights along the trunk at 0.1 m; 0.3 m; 0.5 m; 0.9 m; 1.3 m; 2 m and from 2 m it was measured from meter to meter up to the minimum diameter of 5 cm. The diameter of the logs was measured with a bevel gauge, and lengths with a measuring tape. The individual volume of the trees was determined by Smalian’s formula, according to ^{Machado & Figueiredo Filho (2006}). Diameter measurement was used to obtain individual volumes for posterior comparison with the volume predicted by the equations. A statistical summary of the main dendrometric variables of the 558 sampled trees is described in Table 1.
Statistics | DBH (cm) | h (m) | v (m³) |
---|---|---|---|
Minimum | 14.9 | 14.1 | 0.1374 |
Mean | 37.7 | 29.2 | 1.7473 |
Maximum | 65.0 | 41.5 | 5.4871 |
Standard deviation | 9.1 | 5.0 | 0.9554 |
Coefficient of variation | 24.2 | 17.0 | 54.7 |
Subsequently, five diametric classes were established through previous analysis by the empirical criterion for tree representation.
2.3. Logistic Model
The use of the logistic function to describe tree tapering was addressed by ^{Calegario (2002}) and aims to describe height variations when the diameter also varies (Equation 1). This approach differs from others traditionally used in forest science since (in this case) the diameter variation is the independent variable of the model.
Modified logistic model (^{Carvalho et al., 2014}):
h _{ ij }: height of the i-th tree in the j-th position of the shaft (m); ht _{ i }: total height of the i-th tree (m); r _{ ij }: radius of the i-th tree in the j-th position of the trunk (cm); rap _{ i }: radius of the i-th tree measured at breast height (cm); φ _{ i }: regression of coefficients of fixed and random effects; ε _{ ij }: random error.
The representation of the logistic model is that β are fixed effects, b _{ i } random effects and ε _{ ij } deviations obtained with the prediction of the model regarding the observed variables. Random effects are assumed as being independent for different situations in modeling and the errors within the groups are assumed as being independent for the different scenarios (ε _{ ij } ). The decomposition of model parameters was performed to insert the random effects of diametric class, age and site. The finality is obtaining parsimonious models of simple interpretation and that increase the accuracy of the predictions. The initial estimates of the parameters were obtained through the SelfStart functions, described by ^{Pinheiro & Bates (2000}) through the nlme packet of the Software R.
2.4. Fifth degree polynomial (^{Schöepfer, 1966})
To compare the obtained estimates with the best mixed logistic model, a traditionally addressed model for tapering functions of tree shafts was used (the 5^{th} degree polynomial), as shown in Equation 2.
di: estimated diameter (cm); hi: height along the trunk (m); DBH: diameter at breast height (cm); h: total height (m); β _{ 0 } , β _{ 1 } , β_{ 2 } , β_{ 3 } , β_{ 4 } and β _{ 5 }: parameters of the equation to be estimated; ε _{ ij }: error of the estimate.
2.5. Partial and total volumes
For the volume of the trees, techniques of integration of the base area on the length of the shaft of solids of revolution were used. These techniques are determined through numerical calculation and generate the volume of solid cylindrical shells. This method to obtain the volume of solid shells involves rotating elements of rectangular area parallel to the axis of revolution (axis y), which generates a solid contained between two cylinders with the same center and axis. Thus, the volume of this solid is obtained by the sum of n rectangular elements. As described by ^{Carvalho et al. (2014}), the volume per integration for different radii is obtained with Equation 3 for the modified logistic model.
R _{ min }: radius in the estimate position of minimum radius; R _{ max }: estimate of the maximum radius; R _{ i }: mean radius of the i-th generated cylinder (m); ht _{ i }: total height of the i-th tree (m); r _{ ij }: radius of the i-th tree in the j-th position of the trunk (cm); rap _{ i }: radius of the i-th tree measured at breast height (cm); φ _{ i }: regression coefficients; V: volume of the section between minimum and maximum radius, and consequently the individual volume of the tree.
In order to obtain the volume through the 5^{th} degree polynomial, the sectional area of the tree base was integrated between the lower boundary (h_{1}) and the upper boundary (h_{2}), which is desired to establish as ^{Scolforo (2005}) describes in Equation 4.
v: estimated total or partial volume; h _{ 1 }: lower boundary used in the integration process; h _{ 2 }: upper boundary used in the integration process; w: tapering model as a function of the independent variable d.
In order to assess the accuracy of the estimates obtained with the tapering functions for partial volumes, the tree trunks were divided in three parts (base, middle and apex). The stratification of the shaft was considered, in the base, from 0.1 m to 25% of the total height, in the intermediate portion from 25% to 75% of the total height and, in the apex, from 75% to 95% of the total height (^{Môra, 2015}). The integration was performed using the integrate function associated to the maply function, both implemented in the Software R. For comparison between the best prediction with the 5^{th} degree polynomial, from the logistic model, only the equation that best adjusted between the three random effects analyzed was selected (diametric class, age and site).
2.6. Statistical analysis and model accuracy
To evaluate the proposed models, the Akaike Information Criterion (AIC), the Bayesian Information Criterion (BIC) and the Standard Error of the Estimate (S _{ yx } ) were used, described by ^{Carvalho et al. (2014}). The smaller their respective values, the better the equation was considered.
The graphs of the existing residuals were also analyzed in the prediction of the models, by the difference between the observed values, obtained in the field, and the estimated ones through the equations. The quality of the adjustment of the equations was also assessed through standardized residual dispersion (^{Souza, 1998}). All the equations had their coefficients analyzed through the t-test, with a 5% probability level. The analyses were performed by the Software R, version 3.1, using the nlme package developed by José C. Pinheiro and Douglas Bates (^{R Development Core Team, 2015}).
3. RESULTS AND DISCUSSION
With the use of the theory of mixed models, the fixed and random parameters were estimated by age (Table 2), diametric class (Table 3) and site (Table 4) for the modified logistic model. The significance of the equations was verified in it through the t-test in all the analyzed coefficients, both with fixed and random effects. In the traditional approach, the coefficients were estimated with the 5^{th} degree polynomial (Table 5).
Age | φ1 | φ2 | φ3 | φ4 | |||||||
---|---|---|---|---|---|---|---|---|---|---|---|
Fixed | Random | Fixed | Random | Fixed | Random | Fixed | Random | ||||
11 | 0.8546 | 0.8400 | -0.0667 | -0.0618 | 0.7637 | 0.7372 | 0.1459 | 0.1647 | |||
12 | 0.8141 | -0.0504 | 0.7589 | 0.1545 | |||||||
13 | 0.8562 | -0.0749 | 0.7617 | 0.165 | |||||||
19 | 0.8264 | -0.1047 | 0.8142 | 0.1527 | |||||||
21 | 0.8355 | -0.0696 | 0.762 | 0.1453 | |||||||
23 | 0.8955 | -0.0724 | 0.7386 | 0.1509 | |||||||
24 | 0.8582 | -0.0644 | 0.7513 | 0.142 | |||||||
25 | 0.8721 | -0.0633 | 0.7496 | 0.1484 | |||||||
26 | 0.8425 | -0.0629 | 0.7775 | 0.1402 | |||||||
27 | 0.8606 | -0.0606 | 0.7716 | 0.1313 | |||||||
28 | 0.8727 | -0.0654 | 0.7663 | 0.1351 | |||||||
30 | 0.8823 | -0.0731 | 0.7721 | 0.141 | |||||||
31 | 0.8634 | -0.0712 | 0.7737 | 0.1366 |
Class | φ1 | φ2 | φ3 | φ4 | |||||||
---|---|---|---|---|---|---|---|---|---|---|---|
Fixed | Random | Fixed | Random | Fixed | Random | Fixed | Random | ||||
20 | 0.8647 | 0.8454 | -0.0696 | -0.0608 | 0.7578 | 0.7503 | 0.1492 | 0.1615 | |||
30 | 0.8511 | -0.0647 | 0.7759 | 0.1432 | |||||||
40 | 0.8732 | -0.0694 | 0.7559 | 0.1446 | |||||||
50 | 0.8713 | -0.0758 | 0.7616 | 0.1457 | |||||||
60 | 0.8866 | -0.0958 | 0.7502 | 0.1563 |
Site | φ1 | φ2 | φ3 | φ4 | |||||||
---|---|---|---|---|---|---|---|---|---|---|---|
Fixed | Random | Fixed | Random | Fixed | Random | Fixed | Random | ||||
I | 0.8661 | 0.8810 | -0.0678 | -0.0693 | 0.7642 | 0.7465 | 0.1431 | 0.1525 | |||
II | 0.8578 | -0.0673 | 0.7685 | 0.1421 | |||||||
III | 0.8718 | -0.0663 | 0.7686 | 0.1369 | |||||||
IV | 0.8313 | -0.0465 | 0.7839 | 0.1421 |
Parameter | Estimate | Standard error | t value | Pr (>|t|) |
---|---|---|---|---|
β_{0} | 1.1672 | 0.001547 | 754.7 | <0.0001 |
β_{1} | -3.7173 | 0.041981 | -88.5 | |
β_{2} | 16.7508 | 0.303712 | 55.1 | |
β_{3} | -36.3726 | 0.857805 | -42.4 | |
β_{4} | 34.5127 | 1.035597 | 33.3 | |
β_{5} | -12.3560 | 0.448053 | -27.6 |
In general, in Tables 2, 3 and 4 a small variation is noticed between the parameters of regression. This, in greater magnitude, can contribute to large differences in the prediction process of the settlement assortments.
The statistics of the modified logistic model with the decomposition of the parameters by random effects are shown in Table 6.
Model | AIC | BIC | Syx (%) |
---|---|---|---|
Logistic - Site | -41,169.66 | -41,055.52 | 15.0 |
Logistic - Age | -42,281.13 | -42,166.99 | 14.5 |
Logistic - Diametric Class | -41,225.65 | -41,111.57 | 14.8 |
5^{th} degree polynomial | 63,634 | 63,398 | 1,7 |
In which: AIC is the Akaike Information Criterion; BIC is the Bayesian Information Criterion; and Syx is the residual standard error.
There was no representative distinction between the mixed models for the analyzed statistical criteria. The logistic models demonstrated better adjustments in relation to the 5th degree polynomial for AIC and BIC. For these two statistical criteria used, models with random effects obtained greater accuracy. The logistic model with the random effect age was considered the one with best adjustment in general and, therefore, it was selected for comparison with the 5^{th} degree polynomial. However, for the standard error of the height estimate (hi), one can notice the great improvement that the 5^{th} degree polynomial provided to the estimates. ^{Pires & Calegario (2007}) studied the shaft profile of Eucalyptus through the logistic model and observed standard error values in this range in different categories of diametric classes.
Mixed models were also evaluated through the graphic analysis of the standardized residue for the modeling with random effects age (Figure 1 - A), diametric class (Figure 1 - N) and site (Figure 1 - C) on the hi/ht ratio of the estimated height variable.
The residual distribution found in Figure 1 indicates in general that the residuals in different ages, diametric classes and sites show a similarity between themselves. Through them, it is possible to infer the hypothesis of normality of the data. Thus, because no distribution showed accentuated asymmetry, these can be considered of normal distribution. An important factor is that the data should be distributed next to the zero axis and preferably forming a horizontal square. Again, in this evaluation, these residuals show this aspect. It can be also said that the greater distribution of the residuals in the predictions of the tree shafts, regardless of the random effect, is found in the median region of the trunk. This is noticeable because the greater randomness and values, more distant from the zero axis there are in this portion. Another observation on the model is that it shows a good distribution in the base of the trunk, in which the hi/ht ratio is smaller. Equations with a better residual distribution in the base are preferable, since it is precisely in this point of the tree trunk that logs of higher added value are. ^{Pires & Calegario (2007}) analyzed the shaft profile for logistic models and also found for the tree bases a smaller randomness of residuals. On the other hand, ^{Horle et al. (2010}) found trends of super-estimates for the base of the trunk with the 5^{th} degree polynomial and of the logistic model to describe the profile of Pinus oocarpa. An alternative that could possibly improve these results, and consequently the distribution of residues, would be inserting the effect of their spatial correlation in the modeling.
The residues for the diametric classes in Figure 1 (B) of 20 and 60 cm were those with a more homogeneous distribution. Thus, some outliers were present, especially in diameter classes of 30, 40 and 50 cm. These influence points were remarkable in the initial and median portion of the trunk. Testing the adjustment of the modified logistic model for the tapering of the trunk of Eucalyptus and dividing them in diametric classes, ^{Pires & Calegario (2007}) also verified that this model was superior in the distribution of residues. These authors proved that between linear and non-linear models the logistic one is recommended, since it shows desirable characteristics, such as parameter interpretation, parsimony and data extrapolation. ^{Mendonça et al. (2014}) verified the same behavior of the logistic model to predict settlements assortment of eucalyptus.
The adjustment statistics of the partial volume in the three shaft portions for the two evaluated models are found in Table 7. From this point, only the logistic function with the random effect age was analyzed, since it obtained the best adjustment between the studied effects.
Model | S_{yx-base} (%) | S_{yx-middle} (%) | S_{yx-apex} (%) |
---|---|---|---|
5^{th} Degree Polynomial | 7.9 | 13.7 | 38.7 |
Mixed Modified Logistic - Age | 57.2 | 59.5 | 38.9 |
In which: Syx (%) is the standard error of the estimates.
Through the standard error of the estimates of the partial volumes of the trunk (Table 7), it is noticeable that, regardless of the analyzed portion, the best estimates were produced by the 5^{th} degree polynomial. When verifying this result, the mean standard error of the estimates was analyzed between the portions of the trunk; the 5^{th} degree polynomial obtained 20.1%, and the logistic 51.8%. The 5^{th} degree polynomial was also considered one of the most accurate by ^{Môra (2015}) in the estimates of different portions of the trunk volume for Pinus taeda. However, the errors (S_{yx}) found by this author were less accurate than those observed in this study: 10.8% for the base portion, 20.3% for the intermediate portion, and 51.3% for the final phase of the trunk. As for ^{Souza et al. (2012}), they analyzed the total volume estimates for Pinus taeda and found S_{yx} values smaller than 6% for the equation of the 5^{th} degree polynomial. ^{Téo et al. (2013}) assessed the same model of Pinus elliottii of different ages for the total volume estimates and found errors from 11% to 20.7% (S_{yx}).
Residuals of partial volume estimates from 0.1 meter to 25% of the total height, of the intermediate portion volume from 25% to 75% of the total height, and from the final portion volume from 75% to 95% of the total height are shown in Figure 2.
For the residual produced by volume estimates in the three portions of the tree trunk, the 5^{th} degree polynomial showed the best adjustment again. With graphic analysis, the results obtained in Table 7 are verified, in which the greater accuracy is located in the basal portion of the tree trunks. The predictions obtained with the Mixed Modified Logistic Model had, in the three portions of the shaft, few accurate errors with bias and underestimated, regardless of the diameter of the trees. ^{Môra (2015}) assessed the residuals of the respective partial volumes obtained with the 5^{th} degree polynomial and found, in general, similar results to this study. The volume up to 25% of the total height for trees with DBH < 20 cm was overestimated and the larger ones had constant residual distribution regardless of the DBH. Intermediate and upper portions behaved similarly to the ones observed in Figure 1 (A) and 1 (B). ^{Silva et al. (2011}) verified the accuracy of the diameter estimates along the trunk and the total volume for Pinus caribaea and found that the 5^{th} degree polynomial was again one of the most accurate equations to produce these estimates.