Can an Exponential Function Be Applied to the Asymptotic Density – Size Relationship ? Two New Stand-Density Indices in Mixed-Species Forests

This study presents two stand-density indices (SDIs) based on exponential density decline as a function of quadratic mean diameter for all species combined in mixed-species forests with 22 species mix grouped in four species groups. The exponential-based density–diameter relationship, as well the density index corresponding to the slope or instantaneous mortality rate parameters, was compared with those based on power-law density–diameter relationship. A dataset of 202 fully stocked circular plots at maximum density was used for fitting the models, and a dataset of 122 circular plots was used for validation stand density index for all species combined of mixed-species stands. The dataset for validation was independent of dataset for model development. The first stand-density index showed a density management graphic (DMG) with a variable intercept and common instantaneous mortality rate, and the second index showed a DMG with common intercept and variable mortality rate. Additionally, the value of the initial density of the fitted line was more realistic than those generated by the potential model for all species combined. Moreover, the density management diagrams showed a curvilinear trend based on the maximum stand density index in graphical log–log scale. The DMGs could be interpreted as forest scenarios based on variable initial density and common management objectives or the same density and different management objectives for forest-rotation periods involving all species combined in mixed-species stands. The fitting of exponential and potential equations for species or species groups showed that the density–size relationships in mixed-species forests should be modeled for all species combined because the disaggregation of mixture species represented a weak tendency for each species or species group and the resultant fitted equations were unrealistic.


Introduction
In ecology, density is commonly defined as a number of trees per unit area and competition takes place if the site resources available to each tree are reduced in a given density stage [1].Quantification of site occupancy or stand density is an essential tool for modeling forest stand mortality, growth, and yield [2]. Density is a general concept used to quantify the abundance of trees per unit area in an ecosystem.In forest management, stand density is a term used to describe tree cover or stocking per unit area [3] and this term can be used to relate the tree shape, growth, and mortality.Lately, the stand growth volume has been related with the stand density for making informed management decisions [4], and the stand-density index (SDI) is an important predictor for estimating stand-level biomass [5].The density level is an indicator of forest integrity, particularly because the stand density for a given tree size has a unique limit for a specific species or species group and is independent of other factors as age or site quality [6].The density level is also related to space occupation, which is expressed as the amount of resources used by trees in relation to the maximum resources available on a given site [7].The maximum stand density for a specific species on a given site is an essential element of information for assessing site productivity, modeling and predicting stan dynamics, and designing silvicultural treatments [8].The maximum site occupancy or carrying capacity is a key concept in both ecology and forestry [9].For mixed-species forests, site occupancy is normally defined for all combined species.The regulation of stand density via initial spacing and (or) thinning treatments are, among the controls of stand density, some of the oldest most commonly used methods for achieving forest management goals [4].
In stand-density studies, there are two procedures that have defined the direction of forestry research.The first is based on a research conducted by Reineke [10] for even-aged stands.This approach considers that the number of trees per unit area (N) for full or complete densities varies according to the quadratic mean diameter (QMD) of a stand.Fully stocked stands with smaller QMD have a higher number of trees, whereas fully stocked stands with a larger QMD have relatively fewer trees.The change in N based on a mean diameter, basal area, or QMD follows a regression model with a slope of −1.605 on a log-log scale.The second procedure is based on the "−3/2 power rule of self-thinning" developed by Yoda et al. [11] for pure stands incurring density-dependent mortality.Even though both procedures are algebraically equivalent [12], the Reineke model is based on density-size (number of trees per hectare and the QMD) relationship while Yoda's model on the size-density relationship (mean tree volume and number of trees per hectare).The relationship describes the reciprocal change in biomass or volume and the number of trees per hectare for pure even-aged stands with full density, with the decrease following a model with a slope of −1.5.This relationship is known as the self-thinning rule [6].These two approaches have been the basis for the development of density research in pure and mixed-species forests with different degrees of density and spatial structure [1,[13][14][15][16][17][18].Relative spacing as the average distance between trees divided by the average height of dominant canopy is also used as a measurement of the size-density relationship [2].
The two mean procedures have been a point of discussion in forestry research, with studies both for and against these density approaches.Lonsdale [19] reported that there was no evidence to support the self-thinning rule, because the slope can vary considerably when data is inconsistent.Additionally, the author concluded that the self-thinning line should be defined under controlled forestry conditions, and that fluctuations in the level of available resources at the site should alternate the intercept and consequently the slope in controlled experiments.Hamilton et al. [20] defended the self-thinning rule by noting that the size-density trajectories followed the self-thinning lines in plant populations, but did not necessarily follow a slope of −1.5.Cao et al. [21] showed a density model based on the curvilinear relationship of QMD and density, assuming variable mortality states as density decreases.Ducey and Larson [22] reported that the original form of the SDI defined by Reineke [10] is unrealistic, and that an additive version with the basal area would represent the correct form for different forest-management purposes.Furthermore, Cao and Dean [16] used segmented regression to model the density trajectories for individual stands and QMD over time, with the segmented model characterizing three mortality states.The maximum density or space occupancy can be different between species, and a given absolute value can be represented by a relative density for each species.Therefore, the maximum density limit for all species combined, usually in terms of number of trees or basal area, might be more appropriate for mixed-species forests [3].There is currently a worldwide trend in management of mixed-species stands, and publications on density management have increased considerably [3,8,14,18,[23][24][25]. Maximum stand-density index (SDImax) is an important factor controlling stand dynamics that varies by species and can be used to assess full site occupancy based on species composition [6,26].
The density management diagram (DMD) is a graphical tool for relating stand density, tree size, and stand yield.This represents mean tree volume and stand density.Some of the following relationships have been superimposed in a graph: maximum size-density relationship, imminent competition mortality, crown closure, estimations of diameter and height, and relative density index [27][28][29].The DMDs are used to quickly examine alternative density management regimens and are based upon several ecological and silvicultural concepts [30].
In forest-management planning of pure and mixed-species forest stands in Mexico, power-law density-diameter relationships based on the Reineke model [10] and, occasionally, the self-thinning rule based on Yoda's model [11] have been used to characterize density.Most research has been based on the fitting of size-density or density-size relationships and the construction of DMDs to prescribe thinning [15,17,18,31].Torres-Rojo et al. [32] identified specific details regarding the history of forest management in Mexico and the forest-regulation methods used.In Mexican forests, the SDI based on the first relationship has been used for pure and mixed-species stands based on a reference quadratic mean diameter (QMDR) of 25 cm.However, the power-law density-diameter relationship does not characterize the asymptotic size-density relationship for the different mortality states that occur in a specific stand, and the intercept value of the line fitted to the density axis illustrates the problem of the intercept not being realistic.The hypothesis of this study considers that an exponential based density-diameter relationship can be applied to model the asymptotic maximum density-size relationship in mixed-species forests.The objectives of this study were to develop two modified SDIs based on an exponential relationship between the number of trees per hectare (N) and the QMD and compare the fitting of the equations and SDIs with those based on the density relationship of Reineke's model for all species combined in mixed-species stands under forest management in Durango, Mexico.

Study Area and Data
The dataset used in this study was collected in mixed-species stands of the Forest Management Unit (UMAFOR 1005) "Santiago Papasquiaro y Anexos" in Durango, Mexico.The UMAFOR 1005 is located in the northwestern region of the state of Durango between the coordinates 24 • 30 -25 • 27 N and 105 • 01 -106 • 24 W. The climate is temperate and temperate subhumid with long winters.The average annual temperature is 5 • C in winter months and 27 • C in summer months.The hottest period is from May to June, with a temperatures range of 26-28 • C, and the coldest period lasts from December to February, with temperatures between −6 and 4 • C [33,34].The size of the area with timber production is ~120,014.95ha.
In self-thinning or upper boundary line studies there is a lack of objectivity in determining which data points to include in the fitting of maximum line density [35].This subjectivity in data selection should be avoided with methods such as stochastic frontier regression [36] or quantile regression [26,37].The dataset in this study considered 202 circular plots of 0.10 ha with maximum density, which were selected from 25,017 plots from the forest inventory.In order to avoid the subjectivity in the selection of dataset [35] the experimental data was selected according the following procedure: (1) for all species combined, the density per hectare (N, trees ha −1 ) was plotted against quadratic mean diameter QMD (cm) for visual review; (2) N was ranked by percentile rank for each QMD class (i.e., QMD classes were considered from 13 to 56 cm in 1-cm intervals) in the overall dataset; (3) the 95th percentile was used to select the fully stocked plots of the N-QMD relationship for each QMD class; (4) the resultant N plotted against QMD was used to stablish the density-size relationships for all species combined.This relationship was also fitted by species and species groups.
The information came from mixed-species stands with , and Pse: Prunus serotina Ehrh).The dataset details for the basal area (G, m 2 ), total volume (Vt, m 3 ha −1 ), quadratic mean diameter (QMD, cm), number of trees per hectare (N, trees ha −1 ), SDI developed and estimated with Equation (2) (EE-SDI (β0) ) and SDI for Reineke's model [10] (PE-SDI (β0) , Equation (7)) are shown in Table 1.The dataset using in the fitting process did not evidence density-dependent mortality.Therefore, the maximum stand density was studied for all species combined in mixed-species forests and for species or species group.2) and ( 7)), respectively.
In the 202 plots, 12,684 trees were measured, with P. durangensis having the highest density percentage (24.11%) in the dataset and P. tremuloides the lowest percentage (0.10%).The summary of variables for 22 species is presented in Table 2 for QMD and N. Also, Table 2 shows the number of plots when each species is presented.P. arizonica and P. durangensis were represented in 163 and 161 of the 202 mixed-plots dataset, respectively.

Model Derivation
Assuming that the number of trees per hectare, with an initial density, varies inversely according to the QMD, the rate of change will be inversely proportional to the size.Thus, this relationship can be modeled with a constant rate of mortality expressed by the exponential decline in density (dN/dQMD = −kN) or in the integral form as N 1 = N 0 e −kQMD [38], where N 0 and N 1 are the density at the beginning and end of a given range of QMD, respectively, and k is the instantaneous rate of change or mortality a given range.This approach has been used in ecology and forest studies to evaluate mortality for a given different stage [19,[39][40][41][42][43]; however, this was not studied using the SDI with the equation in integral form.In this study, the asymptotic maximum density-size relationship for all species combined in mixed-species forests was expressed as a density management graphic (DMG) because this relationship is considered key in the construction of DMDs [27,30].Each line on the constructed DMGs represents the corresponding SDI at different percentage for maximum density line, self-thinning line, lower line of constant growth and lower line of free growth or crown closure [28,44,45].Size-density and density-size are fundamental relationships in a DMD [46].
Two new SDIs were developed (QR-SDI-1 and QR-SDI-2) to model the density-size relationship of mixed-species stands when all species were combined.The SDIs are based on the exponential model, and the expression is given by Equation ( 1), with the error term in the additive form.
where N is the number of trees per hectare (trees ha −1 ); QMD is the quadratic mean diameter (cm); β 0 represents the initial density or intercept parameter when QMD is equal zero or close to zero; β 1 represents the instantaneous mortality rate or slope parameter; and ε represents the error associated with the model.The density limit with a reference quadratic mean diameter (QMDR) is used to define the SDImax: SDI = β 0 e (β 1 QMDR) .
Two SDI formulations were developed based on the parameter that was solved from Equation (1).QR-SDI-1 was generated when β 0 is dependent upon the density or SDI (i.e., density-dependent intercept term), and QR-SDI-2 (i.e., density-dependent slope term) was developed for parameter β 1 (Equations ( 2) and (3), respectively): (2) Equation ( 3) can generate a DMG (i.e., maximum density-size line, self-thinning line, lower line of constant growth, and lower line of free growth or crown closure were represented) with a different initial density and an instantaneous mortality rate or a common slope for each density percentage or SDI based on the SDImax, whereas Equation ( 4) can generate a DMG with a common intercept and variable slope.Two expressions estimate the SDI for different combinations of N and QDM at QMDR (i.e., 25 cm, in this case).For a known SDI at QMDR, the number of trees per hectare can be estimated from the corresponding SDI equation (Equations ( 3) and ( 4)) as follows: N (β 1 ) = β0 SDI β0 The proposed equations were compared with the equation generated by Reineke [10] and used in previous studies [9,17,47,48] by means of the nonlinear form represented in Equation (6).The equations for SDIs for a specific QMDR and N when β 0 is dependent upon the density or SDI are given by Equations ( 7) and ( 9), respectively.When β 1 is dependent upon the density or SDI, the expressions are generated as in Equations ( 8) and (10).
where β 0 represents the intercept parameter when QMD is equal to 1; β 1 represents the parameter of the slope; log is the natural logarithm function; and ε represents the error in an additive form.
The disaggregation of mixed-species SDI can be applied with the density proportion of the ith species [15,18].For example, if Equation ( 2) is used, the SDI for a specific species (SDI i ) can be obtained with Equation (11).
where SDI MS is the mixed-species SDI for all species combined; PS i is the ratio or proportion of the ith species and all combined species (i.e., the density, in this case); and SDI MS = ∑ n i N e − β1 (QMD−QMDR) PS i where n is the number of species in mixed-species stand.

Model Fitting and Evaluation
The species were grouped in six species groups: (1) Pinus species group, (2) Quercus species group, (3) other conifers species group, (4) other broadleaves, (5) P. arizonica, P. ayacahuite, P. durangensis, and P. teocote species, and (6) P. durangensis species.The nonlinear least squares method of the statistical environmental R language [49] was used to fit the exponential and potential models for all combined species and species groups.The exponential equation (EE1; Equation ( 1)) was fitted and compared with the potential equation (PE1; Equation ( 6)) developed by Reineke [10].The statistical accuracy of the models was evaluated according to root mean square error (RMSE), the adjusted coefficient of determination (R 2 ), the Akaike information criterion (AIC), the absolute average error (Bias), and coefficient of variation (CV).Additionally, graphical analysis of the curves and fitted lines was considered.The expressions of these statistics can be presented as follows: where N, N and N are the measured, estimated, and average values of the number of trees per hectare (trees ha −1 ), respectively; n is the total number of observations; and p is the number of model parameters.
Fitting equations were validated or tested using a dataset of 122 circular plots (0.10 ha), which were selected from the general dataset, for a quadratic mean diameter close to 25 cm, and independent of the 202 fully stocked plots.Then, the density of each plot was assumed as the theoretical SDI, and the estimated SDI with the corresponding fitted equation was compared to determine the accuracy of the exponential and potential equations.Although this was carried out for all species combined and for each species group, only the information for all species combined is presented.

Results
The estimated parameters for the exponential and potential equations and the fitting statistics are shown in Table 3.The estimate parameters for the exponential and potential equations by species group and fitting statistics are shown in Table 4.The fitting statistics for each species group showed that the size-density relationship was underperforming and some parameters are significantly equal to zero at α = 0.05.The maximum density lines fitted by the two equations to the data for all species combined in mixed-species stands are shown in Figure 1, and Figure 2 shows them on a logarithmic scale (log-log scale).The potential model considers very large intercept values; therefore, the lines of maximum density were projected to a different trend than the data, with these only describing the direction of some points in the dataset.
The maximum density lines fitted by the two equations to the data for all species combined in 262 mixed-species stands are shown in Figure 1, and Figure 2 shows them on a logarithmic scale (log-log

274
(2) and ( 3), two groups of curves associated with SDIs or DMGs were generated to characterize the 275 density for all species combined of mixed-species stands.A DMG with an initial density or variable

276
intercept and instantaneous mortality rate or common slope (Figure 3) and a DMG with a common 277 intercept and variable slope (Figure 4) were constructed.The density line at 100% represented the 278 SDImax (1040), the self-thinning line was defined at 70% density (SDI = 728), the lower limit of the 279 constant growth zone was defined at 40% (SDI = 416), and the lower limit of the free growth zone 280 was defined at 20% (SDI = 260).These results coincided with the definition of growth zones associated

281
with a DMD based on Reineke's model [18] and the Langsaeter theory [50,51].With the exponential equation (EE1, Equation ( 1)), and the estimated SDImax using Equations ( 2) and (3), two groups of curves associated with SDIs or DMGs were generated to characterize the density for all species combined of mixed-species stands.A DMG with an initial density or variable intercept and instantaneous mortality rate or common slope (Figure 3) and a DMG with a common intercept and variable slope (Figure 4) were constructed.The density line at 100% represented the SDImax (1040), the self-thinning line was defined at 70% density (SDI = 728), the lower limit of the constant growth zone was defined at 40% (SDI = 416), and the lower limit of the free growth zone was defined at 20% (SDI = 260).These results coincided with the definition of growth zones associated with a DMD based on Reineke's model [18] and the Langsaeter theory [50,51].
was defined at 20% (SDI = 260).These results coincided with the definition of growth zones associated 281 with a DMD based on Reineke's model [18] and the Langsaeter theory [50,51].

298
The trends of Bias according to SDI class, presented in the form of box-and-whisker plots in Figure 5,

299
showed that for Equation (1), (EE1 and SDI(β0)), the averages of the residuals according to SDI classes 300 were more homogeneous than the other three equations and closest to the zero line.The greatest 301 dispersion was observed in SDI classes of 550 and 750; however, this pattern was similar to the other 302 equations shown.Table 5 shows the validation statistics for the equations fitted to all species combined of mixed-species stands using the exponential and potential equations.The validation dataset revealed the prediction performance of the fitted equations, EE relative to PE.The exponential equation (EE, Equation (1)) produced the most suitable results when Equation (2) was used on the validation data (minimum, average, and maximum QMD values of 24.51 cm, 25.06 cm, and 25.50 cm, respectively) and represented theoretical SDI.Validation statistics based on RMSE, R 2 *, Bias, and CV were better (9.619, 99.378, 1.652, and 1.580, respectively) than those generated by the other implemented SDIs.The trends of Bias according to SDI class, presented in the form of box-and-whisker plots in Figure 5, showed that for Equation ( 1), (EE1 and SDI (β0) ), the averages of the residuals according to SDI classes were more homogeneous than the other three equations and closest to the zero line.The greatest dispersion was observed in SDI classes of 550 and 750; however, this pattern was similar to the other equations shown.2) and (3), respectively); PE-SDI (β0) and PE-SDI (β1) = SDIs based in potential equation (Equations ( 7) and ( 8

317
The SDIs for all species combined and disaggregation of species for 10 randomly selected plots 318 of the fitting dataset is shown in Figure 6.The disaggregation of SDI for species was carried out with 319 Equation ( 11) as a ratio between the number of trees per hectare by species and for all species 320 combined.The number of species belonging to each plot is 6-10 species in mixed-species plots.The 321 P. durangensis and P. arizonica are the main species in those plots.In all cases, the SDIs with the 322 exponential equations (i.e., the SDIs based on a density-dependent parameter and based on a density-

323
independent slope parameter) are greater than those of the potential equation for all species 324 combined and for each species.In some cases, the SDIs when the slope parameter is dependent on 325 density are greater than the corresponding SDI as well as when the intercept parameter is dependent 326 on density (Equations ( 3) and ( 8) for exponential and potential models, respectively).
The SDIs for all species combined and disaggregation of species for 10 randomly selected plots of the fitting dataset is shown in Figure 6.The disaggregation of SDI for species was carried out with Equation ( 11) as a ratio between the number of trees per hectare by species and for all species combined.The number of species belonging to each plot is 6-10 species in mixed-species plots.The P. durangensis and P. arizonica are the main species in those plots.In all cases, the SDIs with the exponential equations (i.e., the SDIs based on a density-dependent parameter and based on a density-independent slope parameter) are greater than those of the potential equation for all species combined and for each species.In some cases, the SDIs when the slope parameter is dependent on density are greater than the corresponding SDI as well as when the intercept parameter is dependent on density (Equations ( 3) and ( 8) for exponential and potential models, respectively).2) and (3), respectively); PE-SDI (β0) and PE-SDI (β1) = SDIs based in potential equation (Equations ( 7) and ( 8), respectively).
For instance, Figure 7 shows the DMGs for the exponential equation with a thinning schedule for 10 mixed-species stands where mixed-species plots were collected.Overall combinations of N and QMD were observed as 70%-100% of SDImax for both DMGs based on a density-dependent parameter (DMG (β0) ) and based on density-dependent slope parameter (DMG (β1) ).For both cases DMG (β0) and DMG (β1) , the thinning schedule was carried out in a constant growth zone (40%-70% of SDImax).The thinning schedule for the mixed-species forest with both DMGs is presented in Table 6.The DMG based on a density-dependent intercept parameter suggests greater removal volume or basal area than the DMG based on a density-independent slope parameter.The average removal percentage of volume or basal area per hectare for the first case was 45.8% while, for the second case, 34.6% was observed for 10 mixed-species stands.6) before thinning schedule and N2 represents the trees per hectare at a specific QMD after thinning schedule with the same QMD for N1 and N2.   2) and (3), respectively) for all species combined; N = number of trees per hectare (trees ha −1 ); BA = basal area (m 2 ha −1 ); V = volume (m 3 ha −1 ); DMG (β0) and DMG (β1) = density management graphics based in exponential equation (Equations ( 2) and (3), respectively).

Discussion
The proposed equation yielded better statistics compared with those generated by Reineke's equation [10], with higher R 2 values and lower RMSE, AIC, and Bias values for the exponential equation compared with those for the potential equation for all species combined.Although the R 2 statistic is used for evaluating linear equations, this was considered as a reference in the fitting and validation process.The exponential equation (EE; Equation (1)) had an R 2 value of 0.9611, which was higher than those reported for Reineke's and Yoda's models on density studies of pure and mixed-species stands [9,17,18].The intercept-estimated parameters or initial densities of the fitted lines for EE and expressed as number of trees per hectare (e.g., β0 = 3120.32)were more realistic than those fitted by PE (e.g., β0 = 44,420.53).EE showed realistic values for the trend of all combined species in mixed-species stand data based on the magnitude of the dataset used to select the plots with maximum density for all species combined, which was attributed to the initial density estimated using the exponential model being very efficient.The initial density using the exponential model assumed a QMD = 0 cm.This value represented the number of trees per hectare when they had not reached a measurable diameter at breast height (dbh) (i.e., when the trees reached 1.3 m, the dbh should be equal 0 cm).The fitted equations for all combined species exposed that the maximum density-size relationship in mixed-species forests should be modeled for all species combined.The density-size relationship when the mixture species are disaggregated by species group or species does not represent the maximum density line (Table 4), because in most of the species or species groups the fitting statistics were poor (i.e., high values of RMSE, Bias and AIC, and low values of R 2 ).Also, in some cases the estimated parameters were significantly equal to zero at α = 0.05.The Pinus durangensis species group (species group 6) showed the best results for species groups for both EE and PE equations, but the fitting statistics were poor.The exponential equation showed greater disaggregated SDI for all species combined or for each species group than the potential equation (Figure 6).Also, in most cases the SDIs calculated with the exponential equation were greater than those estimated by the potential equation.
The slope parameter or instantaneous mortality rate for the exponential equation for all combined species in mixed-species forests ( β1 = −0.0532)was more pronounced as compared to the potential equation based on Reineke's model ( β1 = −1.2698),and this is different to the theoretical value of −1.605 in log-log scale proposed by Reineke [10].In several studies, the slope parameter value has been reported to be around −1.605, with ordinary least square method or stochastic frontier regression for pure and mixed-species forests [2,9,15,17,18,23,52], with and without thinning treatments.This value depends on the species and the species mixed in a given density condition and sometimes takes values around −2.000 [53].The proposed equation for all species combined can be interpreted in the same way as the self-thinning rule [10,11] for the parameter representing the instantaneous mortality rate of the fitted line.The constant mortality rate for all species combined in mixed-species stands suggests a theoretical value of −0.0532 according to the approach using the proposed equation and based on the exponential decrease of the initial density [38] along with changes in QMD classes.The fitting of the maximum density lines was comparable to the lines fitted through quantile or frontier regressions [17,18,31,37] because the data used for the fitting represented the maximum density found in the studied fully stocked mixed-species stands for all species combined (e.g., mean volume or basal area was not used because there was no evidence of density-dependent mortality).The dataset selection was carried out according to the 95th percentile and this guaranteed that the maximum density line was estimated with objectivity [35] for all species combined in mixed-species forests.
The predictive trends of the exponential model were better fitted to the trajectories of the data; the values of the parameter representing the initial density as QMD approaches zero were more realistic and followed a reverse-j form on the normal scale and a curvilinear relationship on the log-log scale, (i.e., in graphical form), that characterized the different episodes of competition in the experimental data (Figure 2).This agrees with the concept proposed by Fang and Bailey [54] for height-diameter models.The average height is equal to 1.3 m when the dbh is equal to 0, that is the QMD should be equal 0 when the tree reaches the height of 1.3 m, and, in this case, the exponential equation represents the initial density at that point.This agrees with results reported by Cao et al. [21], for a specific relationship between QMD and N on the log-log scale, with a linear trajectory for the self-thinning curve for high and curvilinear density levels when densities decreased.
Additionally, this agrees with the procedure used to simulate net forest growth, when curvilinear lines for size-density relationships were assumed [55].A similar approach, but with a segmented model, characterized three episodes of mortality.The first segment represented the initial establishment conditions of the stand, when mortality was not expected, and the other two segments represented two mortality patterns in the data trajectory [16].Graphical analysis confirmed the fitting and prediction statistics of the exponential model, which agreed with the statistics shown in Table 5, and the maximum density line follows a curvilinear form on log-log graphic scale.
The maximum density lines of the DMGs represented the curvilinear tendencies of the SDImax for the different states of mortality or competition associated with the experimental data [16] for all species combined.The DMGs only represented all species combined in mixed-species forests because the maximum density lines for species groups were not realistic (Table 4).The DMGs can be used to define forest scenario strategies or thinning schedules.Thus, a DMG based on a different intercept and common slope can be potentially interpreted as variable initial densities and a similar management objective or, at least, approximately similar during forest rotation (Figure 3).In contrast, a DMG based on a common intercept and different slope (Figure 4) can be interpreted as common initial densities and a different thinning schedule for mixed-species stands.The DMGs based on a density-dependent intercept parameter and a density-dependent slope parameter showed different thinning schedules.The first suggests greater average removal of trees per hectare, volume or basal area than the second (Figure 7 and Table 6).These differences can be associated with the density-dependent parameter computed in each DMG.An important difference in the SDIs developed herein as compared with that proposed by Reineke and the one shown in Equation (7) for the slope parameter is that the proposed SDIs preserve a curvilinear tendency on the log-log scale of maximum density (i.e., in graphical form (Figure 2)), which enables them to integrate aboveground and underground competition, as well as the environmental conditions present in a given site, for a group of species [6] or for all species combined in mixed-species forests as in this case of study.
The parameter for instantaneous rate of change or mortality shown in the integral equation represents the presence of intraspecific competition between trees, referred to as self-tolerance according to a tree-tolerance analogy [56].Tolerant species more effectively use low-intensity light and other more efficient resources relative to intolerant species, allowing them to survive longer in mixed-species stands [6,57].The allometry of the size-density relationship of trees growing under the self-thinning line is particularly informative in regard to eco-physiological aspects and economic timber production, revealing the critical demand for resources in a given growth space [58,59], which is very complex in mixed-species stands.The maximum density can occur after crown closure is full or complete but this is rarely observed because the canopy of a specific stand presents empty spaces and sometimes the size of these spaces is greater than tree growth [18].The condition of density has fluctuations around a level of equilibrium called normal density [60].This condition was represented for all species combined because the fitting statistics in the fitting process were not realistic by species or species groups (Table 4) and the combination of all species represents the site occupancy or carrying capacity [9,61].
The review by del Río et al. [3] described the characterization of the structure, dynamics, and productivity of mixed-species forests and indicated that measurements of absolute density used in pure stands can be used directly in mixed-species stands and should reflect the density patterns of mixed-species stands for all combined species.The SDIs developed in the present study explain the maximum occupation of a given site according to the mixture of species and the density levels at which competition mortality occurs for the species mixed [7,9], and they characterized the SDImax for mixed-species forests (e.g., in this study no evidence of density-dependent mortality was presented and the maximum density-size relationship for species or species groups were not realistic).The maximum density line for all species combined represents the maximum growing space or growing area, the species mixture occupying this growing in both horizontal and vertical directions [3], and the inter-specific and intra-specific competition occurring at the same time [59].The disaggregated SDI by species in each plot or stand can be used in a thinning schedule to regulate the species composition in forest rotation.Natural mortality or self-thinning can be caused by increases in tree size and decreases in self-tolerance, leading to the accumulation of area or gaps between tree crowns [6].Mixing species in a stand can increase the maximum stand density as compared with pure stands under similar site and age conditions [8].To characterize the SDI per species in a given stand, the disaggregation approach of the SDImax can be used for the proportions of species present [15,18] (e.g., Equation ( 11)).The exponential (EE; Equation (1)) and the corresponding SDI when β 0 depends on density (Equation (2)) showed important advantages according to statistical and graphical comparisons.Furthermore, representation of the DMG on the graphical log-log scale generated a maximum density line in curvilinear form that followed the trend of all species combined in the mixed-stand data plots.The developed DMGs can provide resource managers with an objective method of determining a density control schedule by management objectives, and decision-making throughout thinning treatments [28].The theoretical thinning schedule showed that the DMG based on a density-dependent intercept parameter suggests greater average removal percentage in trees per hectare, volume or basal area than the corresponding DMG based on a density-dependent slope parameter (Figure 7 and Table 6).These DMGs with corresponding SDIs can be used in the decision-making process in forestry of mixed-species forests.The SDI lines of each DMG can be used to construct DMDs based on fundamental assumptions about the influence of density on competition, site occupancy and self-thinning [29].Also, relevant ecological and allometric relationships such as yield-density effect and site index [46], and volume or basal area [28,30,45,62] can be included in the DMD.The DMGs constructed with SDImax can be used for forest planning, particularly in determining the optimal timing and intensity of thinning [26].The self-thinning line characterized in each DMG is based on SDImax for all species combined in mixed-species forests and this line represents the ecological limit on the number of trees than can be supported in stands of a certain average size [63].

Conclusions
In this study, the relationship between the number of trees per hectare and quadratic mean diameter was assumed to be exponential for all species combined in mixed-species stands.The fitted maximum density lines for species or species group did not represent the realistic size-density relationship for both exponential and potential equations.Therefore, the maximum density line for mixed-species forest should be considered for all species combined.Additionally, the fitting of the maximum density line for all species combined represented a curvilinear shape on the graphical log-log scale, with a reasonable intercept to the trajectory of the data used along with the exponential decrease in density.The proposed equations based on the exponential equation exhibited better statistical precision than the potential equation of Reineke's model [10] in both all species combined and species groups in mixed-species forests.The SDIs generated two DMGs, with the first showing variable initial density and a common mortality rate, and the second showing a common intercept and a variable slope when the density-dependent intercept or the instantaneous mortality rate parameters were assumed, respectively.The exponential equation and SDIs dependent upon the intercept parameter or initial density yielded the best statistics after model development and validation.The DMGs can be used for thinning prescriptions under two approaches: stands with different initial densities and similar management objectives, or stands with equal initial densities and different management objectives for all species combined in the studied mixed-species stands.

Figure 1 .
Figure 1.Maximum density lines fitted for all species combines of mixed-species plots.EE = the

Figure 3 .
Figure 3. Density management graphic (DMG) with variable initial density and instantaneous

QMDFigure 7 .
Figure 7. Graphical thinning schedule for ten mixed-species stands with DMGs based in exponential equation: SDI dependent on intercept parameter (EE-SDI (β0) ; (a)) and on slope parameter (EE-SDI (β1) ; (b)).N1 represents the trees per hectare at a specific QMD (e.g., Table6) before thinning schedule and N2 represents the trees per hectare at a specific QMD after thinning schedule with the same QMD for N1 and N2.

Table 3 .
Parameter estimates and fitting statistics of the exponential (EE) and potential (PE) equations.

Table 4 .
Parameter estimates and fitting statistics of the exponential (EE) and potential (PE) equations by species groups.Quercus species group, 4 = other broadleaved species group, 5 = Pa, Pay, Pd, and Pt combined species group, 6 = Pd species group; SE = standard error of the estimated parameter; t = Student's t-test; Pr >|t| = probability associated with Student's t distribution; RMSE = root mean square error; R 2 = adjusted coefficient of determination; AIC = Akaike information criterion; Bias = absolute average error.

Table 5 .
Validation statistics of the exponential and potential equations for each SDI.

Table 5 .
Validation statistics of the exponential and potential equations for each SDI.

Table 6 .
Thinning schedule for mixed-species stands with DMGs for the exponential equation.

Table 6 .
Thinning schedule for mixed-species stands with DMGs for the exponential equation.