Individual Tree Diameter and Height Growth Models for 30 Tree Species in Mixed-Species and Uneven-Aged Forests of Mexico

Article Individual Tree Diameter and Height Growth Models for 30 Tree Species in Mixed-Species and Uneven-Aged Forests of Mexico Jaime Briseño-Reyes 1, JoséJavier Corral-Rivas 2,* , Raúl Solis-Moreno 2, Jaime Roberto Padilla-Martínez 2, Daniel JoséVega-Nieva 2, Pablito Marcelo López-Serrano 3, Benedicto Vargas-Larreta 4, Ulises Diéguez-Aranda 5, Gerónimo Quiñonez-Barraza 6and Carlos Antonio López-Sánchez 7 1 Programa Institucional de Doctorado en Ciencias Agropecuarias y Forestales, Universidad Ju á rez del Estado de Durango, Constitución 404 sur Zona Centro, Durango 34000, Mexico; [email protected] 2Facultad de Ciencias Forestales, Universidad Juárez del Estado de Durango, Río Papaloapan y Blvd. Durango S/N Col. Valle del Sur, Durango 34120, Mexico; [email protected] (R.S.-M.); [email protected] (J.R.P.-M.); [email protected] (D.J.V.-N.) 3 Instituto de Silvicultura e Industria de la Madera, Universidad Ju á rez del Estado de Durango, Boulevard del Guadiana 501, Ciudad Universitaria, Torre de Investigación, Durango 34120, Mexico; [email protected] 4 Tecnol ó gico Nacional de M é xico/Instituto Tecnol ó gico de El Salto. Tecnol ó gico 101, Col. La Forestal, El Salto, Durango 34942, Mexico; [email protected] 5 Unidad de Gesti ó n Ambiental y Forestal Sostenible, Departamento de Ingenier í a Agroforestal, Universidad de Santiago de Compostela. Escuela Politécnica Superior, C/Benigno Ledo, Campus Terra, 27002 Lugo, Spain; [email protected] 6 Centro de Investigaci ó n Regional Norte Centro, Campo Experimental Valle del Guadiana, Instituto Nacional de Investigaciones Forestales, Agr í colas y Pecuarias, Carretera Durango-Mezquital Km 4.5, Durango 34170, Mexico; [email protected] 7 Departamento de Biolog í a de Organismos y Sistemas, Escuela Polit é cnica de Mieres, Universidad de Oviedo, C\Gonzalo Gutiérrez de Quirós s/n, 33600 Mieres (Asturias), Spain; [email protected] *Correspondence: [email protected]; Tel.: +52-618-122-9423 Received: 4 March 2020; Accepted: 6 April 2020; Published: 9 April 2020   Abstract: Lack of knowledge of individual tree growth in species-rich, mixed forest ecosystems impedes their sustainable management. In this study, species-specific models for predicting individual diameter at breast height (dbh) and total tree height (h) growth were developed for 30 tree species growing in mixed and uneven-aged forest stands in Durango, Mexico. Growth models were also developed for all pine, all oaks, and all other species of the genus Arbutus (strawberry trees). A database of 55,158 trees with remeasurements of dbh and h of a 5-year growth period was used to develop the models. The data were collected from 217 stem-mapped plots located in the Sierra Madre Occidental (Mexico). Weighted regression was used to remove heteroscedasticity from the species-specific dbh and h growth models using a power function of the tree size independent variables. The final models developed in the present study to predict dbh and total tree height growth included size variables, site factors, and competition variables in their formulation. The developed models fitted the data well and explained between 98 and 99% and of the observed variation of dbh, and between 77 and 98% of the observed variation of total tree height for the studied species and groups of species. The developed models can be used for estimating the individual dbh and h growth for the analyzed species and can be integrated in decision support tools for management planning in these mixed forest ecosystems. Keywords: dbh; total tree height; individual tree growth; species-rich forests Forests 2020,11, 429; doi:10.3390/f11040429 www.mdpi.com/journal/forests Forests 2020,11, 429 2 of 16 1. Introduction Accurate estimates of both current resource levels and the expected resource changes from implementing selected management alternatives are needed to develop tools to support forest management informed decision-making (e.g., [ 1 ]). Forest growth models describe the change in size of a tree or another tree attribute or stand-level variable with age. They are used to analyze and estimate different relationships in forest stand development such as species composition, site characteristics, tree competitive status, and silvicultural management (e.g., [ 2 – 4 ]). Such models, together with stand regeneration, harvesting, and mortality models, are important tools in long-term forest management systems because they can simulate stand development and production under various management alternatives (e.g., [3,5–7]). Empirical forest tree growth models predict growth as a function of tree attributes (such as tree size, crown ratio, vigor and tree competition) and plot or stand level variables (such as age site index and stand density) [ 8 – 12 ]. In particular, individual-tree diameter and height growth models are a fundamental component of forest growth and yield frameworks and constitute very important tools for sustainable forest management planning [ 5 , 7 ]. Although many factors such as tree size, neighborhood competition and environmental variables, have been proposed as being important in explaining patterns of tree growth, their relative contribution for different species and ecosystems are still subject to debate (e.g., [ 5 , 13 – 15 ]). Improving our knowledge of the controls influencing tree growth is central to better understand forest community organization and dynamics [ 15 ], and is also important to the development of sustainable management of forest ecosystems [16]. Understanding the factors that should be included into growth models to explain patterns in mixed-species forest with high species richness is critical for modeling ecological dynamics in such stands [ 16 ]. In uneven-aged, mixed-species natural forests, tree-to-tree variation in growth can be considerable within and among species, as well as spatially within a stand (e.g., [ 17 ]). This variability results from the structural and compositional complexity that often characterize these forests (e.g., [13,17–19] ). In this context, an improved understanding of the factors regulating variability in growth rates can lead to better anticipation of individual tree success and community dynamics [ 17 ]. Such understanding can help in developing biometric tools to support forest management decision-making in mixed-species stands, including the development of management strategies based on natural stand dynamics (e.g., [ 18 , 20 , 21 ]). Nevertheless, factors influencing tree growth in structurally complex forests remain poorly understood [ 18 ]. Most previous studies of competition and tree growth have focused on plantations, managed forests, and even-aged, relatively young natural forests (e.g., [ 22 , 23 ]). In contrast, studies analyzing individual tree growth in mixed forest ecosystems are relatively more scarce (e.g., [ 13 , 17 , 18 , 24 ]), but have received increased attention in the last decades due to a worldwide trend in managing such forests to improve the biodiversity and ecosystems services and to assure long-term sustainable forest resources [3,25]. The Sierra Madre Occidental in the Mexican State of Durango constitutes a hotspot of forest biodiversity, being home to about five million hectares of species-rich and uneven-aged forest stands [ 26 ]. It is the first forest reserve in Mexico with an estimate of 422.1 million m 3 of timber volume, and represents 22.4% of the total estimated production at the national level [ 27 , 28 ]. In 2015, Durango produced 28.5% of the total timber volume harvested in Mexico [ 3 ]. These forests are dominated by irregular and multispecies stands [ 29 ], and have been subject to selective harvesting for more than a century to provide a mix of services to local communities [ 30 ]. Although the main objective of forestry in Durango has been timber production, increased emphasis is being given to the concepts of ecological sustainability and biodiversity. Consequently, the development of forest growth models for mixed and uneven-aged stands is imperative to be used as important tools to support the informed decision-making in this type of forest ecosystems. To date in Durango, forest growth studies have been conducted on site quality [ 31 – 33 ], stand density management [ 28 ], dbh growth models [ 3 ], tree competition [ 34 , 35 ], tree volume and biomass equations [ 27 , 36 , 37 ] and tree diameter-height relationships [ 38 ]. However, most of the previous studies dealt with only one or a few tree species and Forests 2020,11, 429 3 of 16 there have been very limited efforts to develop individual-tree dbh or height growth models for most of the commercial tree species in mixed-species forests. Definitely, we are not aware of studies that have developed individual-tree height growth models for the main tree species in Durango using data of stem-mapped and re-measured plots. In this situation, more research and studies become necessary to better understand stand growth and predict the effect of silvicultural practices in these species-rich and uneven-aged forests that are among the richest in the country. Therefore, the objective of this study was to develop new species-specific models for predicting individual-tree dbh and total tree height growth for 30 commercial species growing in the mixed and uneven-aged stands in Durango, Mexico. Some generic growth equations were also developed for all pines, all oaks, and all strawberry species. Such generic group equations can provide useful estimates for a few remaining tree species with still not enough remeasurements of dbh and total tree height in the area of study. 2. Materials and Methods 2.1. Study Area and Data The data used in this study were collected in the Sierra Madre Occidental in the State of Durango, Northwestern Mexico, located between the geographic coordinates 22 ◦ 20 0 and 26 ◦ 46 0 N latitude and from 103 ◦ 46 0 to 107 ◦ 11 0 W longitude. The data are part of a network of stem-mapped re-measured plots, which were established to monitor the stand dynamics of Durango’s forests [ 38 ]. Plot establishment began in 2007 and re-measurements started in 2012, after a five-year interval. Each plot covers an area of 50 × 50 m. The plots are generally distributed systematically, with a variable grid ranging from 3 to 5 km (Figure 1). The tree number, genus, species, dbh outside bark (dbh; cm), total tree height (h; m), and other variables of all trees with dbh ≥ 7.5 cm were recorded. The installation and re-measurement of the plots followed the methodology developed by Corral-Rivas et al. [ 39 ]. Only undisturbed plots without intermediate harvesting between measurements were selected and considered for this study. A total of 217 re-measured plots were available for the development of the individual-tree dbh and total tree height growth models. The research plots included a variety of stand site conditions, tree sizes and competition levels. A description of the main plot-level variables for the re-measured plots is summarized in Table 1. Forests 2020, 11, x FOR PEER REVIEW 3 of 17 of studies that have developed individual-tree height growth models for the main tree species in Durango using data of stem-mapped and re-measured plots. In this situation, more research and studies become necessary to better understand stand growth and predict the effect of silvicultural practices in these species-rich and uneven-aged forests that are among the richest in the country. Therefore, the objective of this study was to develop new species-specific models for predicting individual-tree dbh and total tree height growth for 30 commercial species growing in the mixed and uneven-aged stands in Durango, Mexico. Some generic growth equations were also developed for all pines, all oaks, and all strawberry species. Such generic group equations can provide useful estimates for a few remaining tree species with still not enough remeasurements of dbh and total tree height in the area of study. 2. Materials and Methods 2.1. Study Area and Data The data used in this study were collected in the Sierra Madre Occidental in the State of Durango, Northwestern Mexico, located between the geographic coordinates 22°20′ and 26°46′ N latitude and from 103°46′ to 107°11′ W longitude. The data are part of a network of stem-mapped re-measured plots, which were established to monitor the stand dynamics of Durango’s forests [38]. Plot establishment began in 2007 and re-measurements started in 2012, after a five-year interval. Each plot covers an area of 50 × 50 m. The plots are generally distributed systematically, with a variable grid ranging from 3 to 5 km (Figure 1). The tree number, genus, species, dbh outside bark (dbh; cm), total tree height (h; m), and other variables of all trees with dbh ≥ 7.5 cm were recorded. The installation and re-measurement of the plots followed the methodology developed by Corral-Rivas et al. [39]. Only undisturbed plots without intermediate harvesting between measurements were selected and considered for this study. A total of 217 re-measured plots were available for the development of the individual-tree dbh and total tree height growth models. The research plots included a variety of stand site conditions, tree sizes and competition levels. A description of the main plot-level variables for the re-measured plots is summarized in Table 1. Figure 1. Study area, and location of 217 stem-mapped plots in Durango, Mexico (plot size is of 50 × 50 m). Figure 1. Study area, and location of 217 stem-mapped plots in Durango, Mexico (plot size is of 50 ×50 m). Forests 2020,11, 429 4 of 16 Table 1. Summary statistics of the 217 re-measured plots used in this study. Variable 2007 2012 Mean Max Min SD Mean Max Min SD N623.78 2264.00 120.00 281.30 624.18 2152.00 144.00 275.91 G21.16 53.63 3.15 8.20 23.31 58.52 3.92 8.98 QMD 21.32 34.10 12.40 4.28 22.34 36.00 13.40 4.44 H016.82 29.50 5.20 4.80 18.18 31.60 5.80 5.19 V188.55 604.21 12.04 106.19 221.04 709.20 16.57 124.28 dbh 17.98 72.52 6.75 9.96 19.19 74.58 7.16 10.26 S8 15 1 2 8 16 1 2 PAI 6.50 21.66 0.24 4.40 N=number of trees per hectare, G=basal area (m 2 ha −1 ), QMD =quadratic mean diameter (cm), H 0 =dominant height, the mean height of the 100 largest-dbh trees per hectare (m), V =timber volume (m 3 ha −1 ), dbh =diameter at breast height (cm), S=number of tree species in plot; PAI =periodic annual volume increment (m 3 ha −1 yr −1 ), SD =standard deviation. The re-measured dbh (dbh 2 ) and total tree height (h 2 ) after a 5-year period were obtained from 55,158 sampled trees of 30 species (each species included more than 100 trees for both variables). The considered species were: Abies durangensis Mart í nez, Alnus firmifolia Fernald, Arbutus arizonica Sargent, A. bicolor Gonz á lez Elizondo, A. madrensis Gonz á lez Elizondo, A. tessellata Sørensen, A. xalapensis Kunth, Cupressus lusitanica Mill., Juniperus deppeana Steud, Pinus arizonica Engelm., P. cooperi C.E. Blanco, P. durangensis Mart í nez, P. engelmannii Mart í nez, P. herrerae Mart í nez Mart í nez, P. leiophylla Schiede, P. lumholtzii B.L. Rob., P. strobiformis Engelm., P. teocote Schltdl., Quercus arizonica Sarg., Q. crassifolia Bonpl., Q. durifolia Seemen, Q. emoryi Porter & J.M. Coult., Q. fulva Liebmann, Q. jonesii Trelease, Q. laeta Liebm., Q. mcvaughii Spellenberg, Q. radiata Trelease, Q. rugosa N é e, Q. sideroxyla Bonpl, Q. urbanii Trel. The information for studied tree species is summarized in Table 2. 2.2. Procedures for Developing the Species-Specific Individual-Tree Growth Equations We used the following empirical model form as starting point of the model selection (1). The equation can be also considered a form of the growth equation proposed by Wykoff[ 40 ] for mixed conifer stands in the Northern Rocky Mountains. It predicts the observed dbh and total tree height for a growth period (5-year in our study) as a function of tree initial size (dbh or height), stand level variables and tree competition as predictor variables. y=∝Tβexp(γS+δC)+ε(1) where yis the dependent variable (dbh —in cm— or total tree height —in m—) re-measured in 2012; T are initial tree size variables (dbh or total tree height) measured in 2007; Srepresents the effect of stand variables that commonly include stand density and dominant height; C represents the effect of tree competition; α , β , γ and δ are the regression coefficients of the equation; and ε is the model error term. Forests 2020,11, 429 5 of 16 Table 2. Summary statistics of the tree species and species groups used to develop the diameter and height growth equations. Species n 2007 2012 Mean SD Mean SD dbh1h1dbh1h1dbh2h2AAI-dbh AAI-h dbh2h2AAI-dbh AAI-h Abies durangensis 134 25.16 16.16 15.05 7.90 26.45 17.74 1.29 1.58 15.25 8.20 0.94 1.65 Alnus firmifolia 425 15.05 9.27 8.32 4.01 16.49 10.86 1.44 1.59 8.63 4.40 1.01 1.50 Arbutus arizonica 592 14.84 5.57 8.09 2.16 15.92 6.32 1.08 0.75 8.34 2.38 0.81 0.87 Arbutus bicolor 698 14.91 5.69 8.24 2.16 16.15 6.61 1.23 0.93 8.69 2.57 0.92 1.05 Arbutus madrensis 846 15.60 5.73 8.14 2.09 16.82 6.62 1.22 0.89 8.56 2.50 0.92 1.06 Arbutus tessellata 387 15.60 5.64 7.96 2.17 16.76 6.51 1.16 0.87 8.18 2.56 0.89 1.24 Arbutus xalapensis 947 14.76 5.90 7.71 2.40 15.81 6.78 1.05 0.88 8.03 2.81 0.91 1.02 Cupressus lusitánica 315 23.78 13.16 19.67 7.95 25.09 14.35 1.32 1.19 20.43 8.47 1.23 1.19 Juniperus deppeana 2896 15.52 7.00 9.02 2.87 16.77 7.82 1.26 0.82 9.29 3.02 0.91 0.87 Pinus arizonica 7618 17.99 10.74 9.27 4.88 19.12 11.95 1.13 1.22 9.59 5.18 0.84 1.07 Pinus cooperi 2770 20.49 13.11 10.25 5.67 22.32 14.52 1.82 1.41 10.79 5.74 1.24 1.18 Pinus durangensis 11,760 18.92 13.23 10.80 5.61 20.33 14.68 1.41 1.45 11.25 5.94 1.06 1.32 Pinus engelmannii 751 23.24 14.56 11.77 6.17 24.65 15.66 1.41 1.10 11.99 6.26 1.14 1.19 Pinus herrerae 1766 25.01 16.62 14.09 6.62 26.54 18.43 1.53 1.81 14.49 6.85 1.23 1.65 Pinus leiophylla 2161 21.09 12.89 10.37 5.54 22.51 14.04 1.43 1.15 10.78 5.80 1.02 1.15 Pinus lumholtzii 1458 16.95 10.50 7.68 4.28 17.98 11.54 1.03 1.04 7.77 4.41 0.77 1.14 Pinus strobiformis 2550 17.74 11.54 9.33 4.72 19.29 12.91 1.54 1.37 9.66 4.92 0.99 1.17 Pinus teocote 3953 19.00 11.95 10.36 5.11 20.43 13.14 1.43 1.18 10.68 5.40 1.00 1.25 Quercus arizonica 179 19.15 9.47 7.79 3.94 20.07 10.91 0.92 1.44 7.95 4.37 0.76 1.58 Quercus crassifolia 2186 17.47 8.14 10.33 3.85 18.46 9.07 0.99 0.93 10.65 4.20 0.87 1.04 Quercus durifolia 255 18.50 7.75 11.02 3.70 19.70 8.44 1.21 0.69 11.18 3.84 0.96 0.91 Quercus emoryi 171 14.54 8.28 8.49 2.63 15.67 9.40 1.13 1.12 8.45 3.01 0.75 1.09 Quercus fulva 327 18.76 9.38 8.58 3.92 19.72 10.42 0.96 1.04 8.65 4.20 0.93 1.25 Quercus jonesii 222 16.19 6.76 8.52 3.74 17.24 7.62 1.04 0.86 8.82 3.98 0.76 0.84 Quercus laeta 389 17.67 9.19 9.64 3.88 18.57 10.28 0.90 1.10 9.82 4.28 0.77 1.31 Quercus mcvaughii 870 14.67 5.96 6.93 2.02 15.72 6.69 1.05 0.73 7.01 2.17 0.73 0.73 Quercus radiata 158 15.65 7.78 8.72 4.04 16.90 8.48 1.25 0.70 9.06 4.25 0.83 0.80 Quercus rugosa 896 18.60 8.78 11.29 3.56 19.97 9.58 1.37 0.80 11.75 3.82 1.05 0.95 Quercus sideroxyla 7373 19.30 9.50 12.26 4.33 20.60 10.56 1.30 1.05 12.65 4.65 1.04 1.20 Quercus urbanii 105 18.97 6.11 8.41 1.79 20.38 6.59 1.41 0.48 8.73 1.92 1.03 0.40 All pines 34,885 19.21 12.45 10.51 5.56 20.59 13.79 1.38 1.33 10.90 5.85 1.03 1.25 All oaks 13,228 18.36 8.82 11.27 4.12 19.56 9.80 1.20 0.98 11.62 4.44 0.98 1.13 All strawberries 3470 15.10 5.73 8.02 2.21 16.25 6.60 1.15 0.87 8.37 2.59 0.90 1.04 n=number of trees, Mean =mean variable values, SD =standard deviation, AAI =average annual increment of dbh (cm) and total tree height (m) observed between 2007 and 2012. Initially, a regression procedure was used to select the definitive tree and stand variables for the individual tree dbh and height growth models. The significance level of 5% was used for the selection of candidate predictive variables in the model. These predictive variables were the measurements at the beginning of the selected 5-year growth period. They included initial tree dbh (dbh 1 , cm); initial tree total height (h 1 , m); the sum of the basal area of the trees larger than the subject tree (BAL, proposed by Wykoffet al. [ 41 ]), a modified BAL (BALMOD, proposed by Schröder and Gadow [ 42 ]) which includes the relative spacing index (i.e., BALMOD =(1 − (1 − (BAL/G))/IH, where IH =100/H 0 N 1/2 ); the stand quadratic mean diameter (QMD, cm); stand density (N, number of trees per hectare); the total basal area of the plot (G, m 2 ha −1 ); stand dominant height (H 0 , m); and crown competition factor (CCF) [ 43 ]. H 0 was calculated from the proportion of the 100 thickest trees of pines per hectare. Data from trees with broken or dead tops were not used to calculate H 0. These predictive variables, along with their various combinations were tested separately by each species or species group using the lm and nls functions of R software environment [ 44 ]. In order to avoid problems related to multicollinearity, the value of the variance inflation factor (VIF) was computed for the groups of variables that best predicted the dependent variables. Candidate variables with a VIF larger than 10 were excluded, as suggested by Myers [ 45 ]. Non-significant parameters at 5% were also removed in the modeling process (i.e., they were fixed at 1 in the final model fitting). Therefore, a little bit different models were used for such species, resulting in the following two dbh and height models as the best candidate model formulations for describing both dbh and height growth in general among the models tested: dbh2=a0dbha1 1exp(a2H0+a3BAL)+ε(2) Forests 2020,11, 429 6 of 16 h2=b0hb1 1exp(b2H0+b3BAL)+ε(3) where dbh 2 is the diameter at breast height after the selected 5-year growth period (cm), dbh 1 is the initial tree diameter (cm), H 0 is the stand dominant height (m), BAL is the sum of the basal area of the trees larger than the subject tree proposed by Wykoffet al. [ 41 ], h 2 is the total tree height after a 5-year growth period (m), h 1 is the initial total tree height (m), a 0 –b 3 are the regression coefficients, and ε is the model error term. The models were first fitted under the assumption that the errors within equation were independent and identically distributed [ 46 ]. As only two measurements were available for each tree in the dataset, autocorrelation among measurements within a tree was not considered in the present study. However, residual analysis revealed the presence of heteroscedasticity in most species-specific dbh and total tree height growth equations. To deal with this problem weighted regression was used, modelling the error variances as a power function of the independent size variables (e.g., [ 46 , 47 ]). For the diameter growth model, the error variance was modelled with a power function of dbh 1 as independent variable (i.e., σ2 i=dbhk 1 ). Using the method suggested by Park [ 48 ], squared residuals ( ˆ ε2 i ) were regressed against dbh 1 to obtain an estimate of kas follows: ˆ ε2 i = γxdbhk 1 , this expression was linearized by the use of the natural logarithm: ln ˆ ε2 i =ln γ+kln dbh1 and parameters were then estimated by linear least squares regression; the species-specific kvalue was subsequently included in the weighting factor ( wi= 1 /dbhk i1 ). The same was done for the species-specific height growth models (3), but using h 1 as the independent variable. 2.3. Model Comparison Evaluation of the performance of the dbh and height growth equations was based on numerical and graphical analysis of the residuals. Two goodness-of-fit statistics were calculated: the root mean square error (RMSE) and the coefficient of determination (R 2 ). These statistics were computed on the basis of the unweighted residuals because these are the values of most interest by users of the equations. The expressions of these statistics are: RMSE =v u tPn i=1yi−ˆyi2 n−p(4) R2=1−Pn i=1(yi−ˆ yi)2 Pn i=1(yi−y)2(5) where yi , ˆ yi and y are the observed, predicted, and mean values of the dependent variable, respectively; and nand pare the total number of observations and the number of parameters used to fit the models, respectively. Plots of studentized residuals against predicted values of the dependent variables were used for visual evaluation of the models [ 49 , 50 ]. These graphs are very important for indicating areas for which the equations provide especially poor or good predictions. 3. Results 3.1. Diameter Growth Equations The parameter estimates and the goodness-of-fit statistics for predicting dbh after a 5-year growth period with model 2 for each analyzed tree species are presented in Table 3. Most of the parameter estimates were significant at α =0.05. All non-significant parameters were dropped from model 2. The initial tree diameter was the most significant and important predictor of dbh 2 . In all cases the coefficients a 0 and a 1 were positive, and a 1 was less than 1.0 in all cases, indicating that dbh 1 limits diameter growth by itself. The term H 0 was positive in all species except for Pinus engelmannii and Forests 2020,11, 429 7 of 16 Quercus arizonica (see a 2 in Table 3), resulting in increasing diameter with increasing stand dominant height. On the other hand, coefficient a 3 multiplying the BAL competition index was negative in all species except for Pinus engelmannii and P. herrerae, indicating a decrease in diameter growth as individual tree competition increases. Table 3. Parameter estimates and goodness-of-fit statistics of the individual tree dbh growth equations developed. Species nEstimated Parameters RMSE (cm) R2 a0a1a2a3 Abies durangensis 134 1.2745 0.9422 ns −0.0013 0.88 0.997 Alnus firmifolia 425 1.3368 0.9453 ns −0.002 0.97 0.987 Arbutus arizonica 592 1.2053 0.9589 ns ns 0.79 0.991 Arbutus bicolor 698 1.2017 0.9619 ns −0.001 0.83 0.991 Arbutus madrensis 846 1.1737 0.9658 0.0013 ns 0.85 0.99 Arbutus tessellata 387 1.2073 0.9626 ns ns 0.87 0.989 Arbutus xalapensis 947 1.1680 0.9721 ns ns 0.87 0.988 Cupressus lusitanica 315 1.1516 0.9735 ns −0.0007 0.98 0.998 Juniperus deppeana 2896 1.2959 0.9421 0.0008 −0.0017 0.87 0.991 Pinus arizonica 7618 1.1711 0.9618 0.0018 −0.0007 0.78 0.993 Pinus cooperi 2770 1.3175 0.9251 0.0045 −0.0027 1.11 0.99 Pinus durangensis 11,760 1.2015 0.9567 0.002 −0.001 0.96 0.993 Pinus engelmannii 751 1.2143 0.9689 −0.0026 0.001 1.12 0.991 Pinus herrerae 1766 1.1624 0.9729 ns 0.0005 1.16 0.994 Pinus leiophylla 2161 1.222 0.9509 0.0021 −0.0014 0.95 0.992 Pinus lumholtzii 1458 1.2734 0.9320 0.0025 −0.0019 0.77 0.99 Pinus strobiformis 2550 1.2577 0.9402 0.0027 −0.0009 0.93 0.991 Pinus teocote 3953 1.2358 0.9431 0.0031 −0.0012 0.94 0.992 Quercus arizonica 179 1.1084 0.9931 −0.003 0.0014 0.75 0.991 Quercus crassifolia 2186 1.1858 0.9598 0.0012 −0.001 0.83 0.994 Quercus durifolia 255 1.3001 0.9386 ns −0.0017 0.94 0.993 Quercus emoryi 171 1.5079 0.882 0.0048 −0.0066 0.62 0.995 Quercus fulva 327 1.3912 0.9165 0.0013 −0.0037 0.91 0.989 Quercus jonesii 222 1.1978 0.9526 ns −0.0026 0.71 0.994 Quercus laeta 389 1.2017 0.9593 ns −0.0011 0.75 0.994 Quercus mcvaughii 870 1.2704 0.9333 0.0023 −0.0012 0.73 0.989 Quercus radiata 158 1.2433 0.9278 0.0069 −0.0026 0.74 0.993 Quercus rugosa 896 1.2257 0.954 0.0015 −0.0013 0.96 0.993 Quercus sideroxyla 7373 1.258 0.9438 0.0022 −0.0018 0.97 0.994 Quercus urbanii 105 1.2524 0.9547 ns ns 0.99 0.987 All pines 34,885 1.2085 0.9516 0.0025 −0.0010 0.96 0.992 All oaks 13,228 1.2328 0.9495 0.0018 −0.0015 0.93 0.994 All strawberries 3470 1.1875 0.9659 ns −0.0004 0.85 0.990 n=number of trees used in model fitting, a 1 –a 3 =parameter estimates, ns =these parameters were not significant at p<0.05, and were therefore fixed at 1 in the final model fitting, resulting in a little different model for these species, RMSE =root mean square error, R2=coefficient of determination. The developed equations fitted the data well and generally explained between 98 and 99% of the observed variation of dbh 2 for all studied tree species. Average RMSE values for pines, other conifers (Abies durangensis, Cupressus lusitanica, and Juniperus deppeana), oaks, and other broadleaf species were 0.93 cm ± 0.17 cm, 0.83 cm ± 0.12 cm, 0.84 cm ± 0.03 cm, and 0.92 cm ± 0.06 cm, respectively. The RMSE values estimated for all pines, all oaks, and all strawberries species groups were 0.96, 0.93 and 0.85 cm, respectively. The graphs of studentized residuals against predicted dbh values showed no evidence of heterogeneous variance over the full range of predicted values for all species after fitting the models with weighted regression. Forests 2020,11, 429 8 of 16 In accordance with Huang et al. [ 50 ], plots of RMSE were obtained across tree sizes to evaluate the accuracy and precision of diameter predictions of model 2 for the developed equations. Figure 2 shows the plots for the three generic equations developed in this paper (species groups: all pines, all oaks, and all strawberries). It can be inferred that the final equations performed quite well for all sampled diameter classes. All strawberries model produced larger prediction errors for large-sized trees in comparison with all oaks and all pines models. Forests 2020, 11, x FOR PEER REVIEW 9 of 17 Quercus sideroxyla 7373 1.258 0.9438 0.0022 −0.0018 0.97 0.994 Quercus urbanii 105 1.2524 0.9547 ns ns 0.99 0.987 All pines 34,885 1.2085 0.9516 0.0025 −0.0010 0.96 0.992 All oaks 13,228 1.2328 0.9495 0.0018 −0.0015 0.93 0.994 All strawberries 3470 1.1875 0.9659 ns −0.0004 0.85 0.990 n = number of trees used in model fitting, a1–a3 = parameter estimates, ns = these parameters were not significant at p < 0.05, and were therefore fixed at 1 in the final model fitting, resulting in a little different model for these species, RMSE = root mean square error, R2 = coefficient of determination. The graphs of studentized residuals against predicted dbh values showed no evidence of heterogeneous variance over the full range of predicted values for all species after fitting the models with weighted regression. In accordance with Huang et al. [50], plots of RMSE were obtained across tree sizes to evaluate the accuracy and precision of diameter predictions of model 2 for the developed equations. Figure 2 shows the plots for the three generic equations developed in this paper (species groups: all pines, all oaks, and all strawberries). It can be inferred that the final equations performed quite well for all sampled diameter classes. All strawberries model produced larger prediction errors for large-sized trees in comparison with all oaks and all pines models. Figure 2. Root mean square error averaged over 5-cm dbh classes for three generic individual diameter growth equations developed in this study. 3.2. Height Growth Equations The parameter estimates and the goodness-of-fit statistics for predicting total tree height after a 5-year growth period with model 3 for each species are presented in Table 4. Most of the parameter estimates were significant at α = 0.05. All non-significant parameters were removed from model 3. The initial total tree height was the most significant and important predictor of height after a 5-year growth period. In all cases, the coefficients b0 and b1 were positive, b1 was less than 1.0 in all cases (except for Quercus urbanii), indicating than h1 limits height growth by itself. The term H0 was positive in all species except for Quercus urbanii and Cupressus lusitanica (see b2 in Table 3), resulting in Figure 2. Root mean square error averaged over 5-cm dbh classes for three generic individual diameter growth equations developed in this study. 3.2. Height Growth Equations The parameter estimates and the goodness-of-fit statistics for predicting total tree height after a 5-year growth period with model 3 for each species are presented in Table 4. Most of the parameter estimates were significant at α =0.05. All non-significant parameters were removed from model 3. The initial total tree height was the most significant and important predictor of height after a 5-year growth period. In all cases, the coefficients b 0 and b 1 were positive, b 1 was less than 1.0 in all cases (except for Quercus urbanii), indicating than h 1 limits height growth by itself. The term H 0 was positive in all species except for Quercus urbanii and Cupressus lusitanica (see b 2 in Table 3), resulting in increasing height growth with increasing stand dominant height. On the other hand, coefficient a 3 multiplying the BAL competition index was negative in all species except for Cupressus lusitanica, indicating a decrease in height growth as individual tree competition increases. The developed equations fitted the data well and generally explained between 77 and 98% of the observed variation of the total tree height for all studied species. Average RMSE values of total tree height of pines, other conifers, oaks, and other broadleaf species were 1.16 m ± 0.25 m, 0.99 m ± 0.29 m, 1.03 m ± 0.12 m, and 1.41 m ± 0.26 m, respectively. The RMSE values estimated for all pine, all oak, and all strawberries species groups were 1.23, 1.10 and 1.03 m, respectively. Forests 2020,11, 429 9 of 16 Table 4. Parameter estimates and goodness-of-fit statistics of the total individual tree height growth equations developed. Species nEstimated Parameters RMSE (m) R2 b0b1b2b3 Abies durangensis 134 1.3983 0.8883 0.0054 ns 1.64 0.960 Alnus firmifolia 425 1.6787 0.8847 ns −0.0019 1.48 0.887 Arbutus arizonica 592 1.2969 0.8920 0.0049 ns 0.87 0.867 Arbutus bicolor 698 1.2665 0.9044 0.0082 −0.0022 1.03 0.840 Arbutus madrensis 846 1.2893 0.9273 ns ns 1.05 0.824 Arbutus tessellata 387 1.5365 0.8387 0.0087 −0.0090 1.22 0.773 Arbutus xalapensis 947 1.2214 0.9441 0.0041 ns 0.99 0.876 Cupressus lusitanica 315 1.2178 0.9759 −0.0018 0.0002 1.10 0.983 Juniperus deppeana 2896 1.3426 0.8885 0.0047 −0.0023 0.87 0.918 Pinus arizonica 7618 1.3660 0.8813 0.0072 −0.0017 1.03 0.961 Pinus cooperi 2770 1.5429 0.8528 0.0049 −0.0021 1.16 0.959 Pinus durangensis 11,760 1.3379 0.9096 0.0034 −0.0006 1.31 0.952 Pinus engelmannii 751 1.3622 0.8881 0.0054 −0.0018 1.19 0.964 Pinus herrerae 1766 1.6105 0.8488 0.0047 −0.0022 1.64 0.943 Pinus leiophylla 2161 1.3041 0.9095 0.0042 −0.0013 1.14 0.962 Pinus lumholtzii 1458 1.3825 0.8874 0.0041 −0.0012 1.15 0.932 Pinus strobiformis 2550 1.3900 0.8924 0.0038 −0.0009 1.16 0.944 Pinus teocote 3953 1.2745 0.9115 0.0051 ns 1.23 0.948 Quercus arizonica 179 1.4913 0.8536 0.0081 −0.0045 1.56 0.872 Quercus crassifolia 2186 1.2162 0.9275 0.0044 ns 1.01 0.943 Quercus durifolia 255 1.2883 0.9023 0.0049 −0.0030 0.90 0.945 Quercus emoryi 171 1.2861 0.9079 0.0078 −0.0039 1.06 0.875 Quercus fulva 327 1.2640 0.8933 0.0074 ns 1.23 0.914 Quercus jonesii 222 1.2470 0.9241 ns ns 0.83 0.956 Quercus laeta 389 1.2588 0.9217 0.0048 ns 1.29 0.909 Quercus mcvaughii 870 1.2981 0.8996 0.0035 ns 0.73 0.887 Quercus radiata 158 1.2473 0.9186 ns −0.0029 0.79 0.966 Quercus rugosa 896 1.2068 0.9290 0.0048 −0.0015 0.93 0.940 Quercus sideroxyla 7373 1.2952 0.9087 0.0046 −0.0014 1.18 0.935 Quercus urbanii 105 1.3727 1.0010 −0.0190 ns 0.37 0.962 All pines 34,885 1.3512 0.8957 0.0050 −0.0010 1.23 0.956 All oaks 13,228 1.2560 0.9156 0.0050 −0.0013 1.10 0.938 All strawberries 3470 1.2800 0.9207 0.0043 −0.0018 1.03 0.843 n=number of trees used in model fitting, b 0 –b 3 =parameter estimates, ns =these parameters were not significant at p<0.05, and were therefore fixed at 1 in the final model fitting, resulting in a little different model for these species, RMSE =root mean square error, R2=coefficient of determination. The graphs of studentized residuals against predicted total height values showed no evidence of heterogeneous variance over the full range of predicted values for all species after fitting the models with weighted regression. Figure 3shows the plots of RMSE obtained across tree height classes for the three generic equations developed in this paper to predict total tree height after a 5-year growth period with model 3. It can be seen that the final models performed quite well for all sampled height classes. The prediction errors for the all strawberries equation can be only evaluated for trees less than 25 m in height due to these tree species occupy always the subdominant vertical layers when mature and are small-sized trees [51]. Forests 2020,11, 429 16 of 16 63. King, D.A.; Davies, S.J.; Nur Supardi, M.N.; Tan, S. Tree growth is related to light interception and wood density in two mixed dipterocarp forests of Malaysia. Funct. Ecol. 2005,19, 445–453. 64. Metz, J.; Seidel, D.; Schall, P.; Scheffer, D.; Schulze, E.D.; Ammer, C. Crown modeling by terrestrial laser scanning as an approach to assess the effect of aboveground intra-and interspecific competition on tree growth. For. Ecol. Manag. 2013,310, 275–288. 65. Tom é , M.; Burkhart, H.E. Distance-dependent competition measures for predicting growth of individual trees. For. Sci. 1989,35, 816–831. 66. Biging, G.S.; Dobbertin, M. Evaluation of competition indices in individual tree growth models. For. Sci. 1995,41, 360–377. 67. Canham, C.D.; LePage, P.T.; Coates, K.D. A neighborhood analysis of canopy tree competition: Effects of shading versus crowding. Can. J. For. Res. 2004,34, 778–787. 68. Coates, K.D.; Canham, C.D.; LePage, P.T. Aboveversus below-ground competitive effects and responses of a guild of temperate tree species. J. Ecol. 2009,97, 118–130. 69. Pretzsch, H.; Bielak, K.; Block, J.; Bruchwald, A.; Dieler, J.; Ehrhart, H.-P.; Kohnle, U.; Nagel, J.; Spellmann, H.; Zasada, M.; et al. Productivity of mixed versus pure stands of oak (Quercus petraea (Matt.) Liebl. and Quercus robur L.) and European beech (Fagus sylvatica L.) along an ecological gradient. Eur. J. For. Res. 2013 ,132, 263–280. 70. Keddy, P.A. Studying Competition. In Competition, 2nd ed.; Usher, M.B., DeAngelis, D.L., Eds.; Springer Science & Business Media: Louisiana, LA, USA, 2001; pp. 1–58. 71. Adame, P.; Hynynen, J.; Canellas, I.; del R í o, M. Individual-tree diameter growth model for rebollo oak (Quercus pyrenaica Willd.) coppices. For. Ecol. Manag. 2008,255, 1011–1022. 72. Pukkala, T.; Lahde, E.; Laiho, O. Growth and yield models for uneven-sized forest stands in Finland. For. Ecol. Manag. 2009,258, 207–216. 73. Á lvarez-Taboada, M.F.; Barrio-Anta, M.; Gorgoso-Varela, J.; Alvarez-Gonz á lez, J.G. Influence of the competition on the section growth in Pinus radiata D. Don. Inv. Agrar.-Sist. Rec. F. 2003,12, 25–35. 74. Corral-Rivas, J.J.; Gonz á lez-Elizondo, M.S.; Lujan-Soto, J.E.; Gadow, K.V. Effects of density and structure on production in the communal forests of the Mexican Sierra Madre Occidental. South. For. J. For. Sci. 2019 ,81, 1–10. [CrossRef] 75. Vargas-Larreta, B.; Castedo-Dorado, F.; Á lvarez-Gonz á lez, J.G.; Barrio-Anta, M.; Cruz-Cobos, F. A generalized height-diameter model with random coefficients for uneven-aged stands in El Salto, Durango (Mexico). Forestry 2009,84, 445–462. [CrossRef] 76. Abetz, P. Zur Standraumregulierung in Mischbeständen und Auswahl von Zukunftsbäumen. AFZ 1974 ,29, 871–873. 77. Monserud, R.A.; Rehfeldt, G.E. Genetic and environmental components of variation of site index in inland Douglas-fir. For. Sci. 1990,36, 1–9. 78. Hasenauer, H.; Burkhart, H.E.; Sterba, H. Variation in potential volume yield of loblolly pine plantations. For. Sci. 1994,40, 162–176. 79. S á nchez-Salguero, R.; Linares, J.C.; Camarero, J.J.; Madriga-Gonzalez, J.; Hevia, A.; Sanchez-Miranda, A.; Ballesteros-Canovas, J.A.; Alfaro-Sanchez, C.; Garcia-Cervignon, A.K.; Bigler, C.; et al. Disentangling the effects of competition and climate on individual tree growth: A retrospective and dynamic approach in Scots pine. For. Ecol. Manag. 2015,358, 12–25. [CrossRef] 80. Bigelow, S.W.; Runkle, J.R.; Oswald, E.M. Competition, Climate, and Size Effects on Radial Growth in an Old-Growth Hemlock Forest. Forests 2020,11, 52. [CrossRef] © 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).