Assessing Differences in Competitive Effects among Tree Species in Central British Columbia, Canada

Research Highlights: We investigated the competitive interactions among three tree species (interior Douglas-fir (Pseudotsuga menziesii var. glauca (Beissn.) Franco), interior spruce (Picea glauca [Moench] Voss × Picea engelmannii Engelm.), and lodgepole pine (Pinus contorta Dougl. Ex Loud. var. latifolia Englem.)) in multi-aged stands in central British Columbia, Canada. Background and Objectives: Understanding competitive interactions among tree species in mixed-species stands is fundamental to supporting silvicultural decision-making in such stands. Using the periodic annual basal area increment for single trees as our dependent variable, we investigated whether neighboring trees competed with subject trees independently of species identity. We also examined the differences in single-tree basal area growth among the three conifer species over time under different levels of competition. Materials and Methods: We developed several spatially explicit, single-tree basal area growth models for interior Douglas-fir, interior spruce, and lodgepole pine using data from 16 plots in two blocks of a long-term study (five measurements over a 21-year period) on the response to pre-commercial thinning. We compared these equations to assess whether intraspecific or interspecific competition predominated. We also examined the differences in basal area growth among the three conifer species over time under different levels of competition. Results: We found asymmetrical relationships between the conifer trees and their neighbors for all species, indicating that the main driver limiting growth in these stands is aboveground competition for light. There was evidence of higher intraspecific competition for small (<10.0 cm DBH) interior Douglas-fir in one block. However, there was no general pattern among larger subject trees with respect to the identity of neighborhood competitive effects and the equivalence of neighbors. We observed a higher level of basal area growth over time for interior Douglas-fir than for lodgepole pine and interior spruce, irrespective of the competition intensity and, not surprisingly, the growth rate declined with increasing competition levels for the three species. Conclusions: Our results provide an understanding of how interior Douglas-fir stands will develop over time and information on species interactions that could help forest managers explore different silvicultural options and their effects on individual tree growth in these complex stands.


Introduction
Mixed-species forest stands are increasingly gaining attention as studies have highlighted their importance in providing higher levels of ecosystem services and functions [1]. Mixed-species stands may also be more productive and have a higher resistance to biotic and abiotic disturbances [2,3]. Interactions between tree species in such forests may be complex and dynamic due to changes to the environment and resource availability. However, understanding the interaction among species in mixed-species stands is fundamental to supporting managerial decision-making. a long-term study (five measurements over a 21-year period) on the response to pre-commer ning in the central interior of BC. These models were used to (i) investigate whether neighbor s compete with a subject tree independently of species identity, and (ii) examine the differen asal area growth among the three conifer tree species present in the study area (interior Dougl nterior spruce, and lodgepole pine (Pinus contorta Dougl. Ex Loud. var. latifolia Englem.) o under different degrees of competition. To address these objectives, we fitted a nonlinear mix t basal area growth model for a subject tree, incorporating tree size (DBH) and the competit t of neighboring trees, irrespective of species (i.e., all species are treated as equival petitors). Subsequently, we partitioned the competitive effects of neighbors into different spec ecies groups to allow us to incorporate potential differences in competitive relationships amo ies.

Study Site and Data
We examined 21 years of growth records from a pre-commercial thinning (spacing) study in cool subzone of the Interior Douglas-fir zone (IDFdk) of the Biogeoclimatic Ecosyst sification (BEC) system in British Columbia, Canada [24]. The study was located in the Kn k Unit of the University of British Columbia's Alex Fraser Research Forest (52°3'N, 121°52'W Knife Creek Unit covers an area of approximately 3500 ha [25], situated on the Fraser plate roximately 27 km east of Williams Lake, in the south-central part of BC ( Figure 1). The Fra eau has an average elevation of about 1000 m, a mean daily temperature of 4.2 °C, and an aver ual precipitation of 450 mm [24]. During the summers of 1989 and 1990, three blocks (B, C, and D), of approximately 40 hecta , were selected for this experiment. Block B is quite dry, Block C is more mesic, and Block D oistest of the sites and is transitional to the Sub-Boreal Spruce (SBS) BEC zone. The stands w  During the summers of 1989 and 1990, three blocks (B, C, and D), of approximately 40 hectares each, were selected for this experiment. Block B is quite dry, Block C is more mesic, and Block D is the moistest of the sites and is transitional to the Sub-Boreal Spruce (SBS) BEC zone. The stands were logged during the 1950s and 1960s to a diameter limit of approximately 25.4 cm, although a few of the larger trees were not harvested. Block B was dominated by interior Douglas-fir (>80% by basal area), with the remaining basal area consisting of lodgepole pine and a small percentage (<0.5%) of trembling aspen (Populus tremuloides Michx). Block C was also dominated by Douglas-fir (>60% basal area). Lodgepole pine varied among the plots (0%-26.5% basal area) and there was a small interior spruce component (less than 5% by basal area). The remaining basal area was comprised of white birch (Betula papyrifera Marsh.) and trembling aspen. Block D was dominated by mixtures of Douglas-fir, spruce, and lodgepole pine (together accounting for more than 80% of the basal area). The remaining basal area was comprised of aspen and white birch. For this analysis, plots from Blocks C and D (16 in total) were used because a sufficient species mixture was present to enable the assessment of interspecific competition in terms of the basal area growth of the three conifer species (hereafter referred to as interior Douglas-fir (FD), interior spruce (SX), and lodgepole pine (PL)) ( Table 1). Table 1. Summary statistics (number of trees (n), and mean, minimum (Min), and maximum (Max) diameter at breast height (DBH)) by block, measurement year, and tree species or species group.

Year
Tree Three pre-commercial thinning treatments were applied during the fall and winter of 1990-1991 to three portions (quarters) of each block and the fourth portion was used as a control. Two 0.05 ha permanent plots were established in dense areas of each of the four portions of each block prior to the treatments, resulting in a total of 24 plots (3 sites × 4 treatments × 2 plots). The thinning treatments employed were the (1) standard approach, where areas were pre-commercially thinned according to the 1990 standards of the BC Ministry of Forests; (2) 3 m clumped approach, where each retained clump included 3-9 trees of the same height class within a 3 m radius circle and the distance between each clump was 3 m; (3) 5 m clumped approach, where each retained clump included 3-9 trees of the same height class within a 3 m radius circle and the distance between each clump was 5 m; and (4) a control, where no thinning treatment was applied [26].
During the spring of 1993, all living trees taller than 1.3 m within the confines of the plots, and those trees greater than 10 cm diameter at breast height (DBH) within a 5 m distance outside the plot boundary, were permanently tagged and measured. Measurements made on tagged trees included the tree location, species, DBH, total tree height, heights to the base of the live crown in each of the four quarters, crown diameter, and tree vigor. Four subsequent measurements (1997, 2004, 2008, and 2013), representing four growth periods, were made. The first four measurements were made prior to the onset of the growing season noted; the 2013 measurements were made following that growing season.
Block C had a higher number of FD than Block D, with a slightly lower mean DBH (Table 1). Mean DBH values for living FD, SX, and broadleaf trees (aspen and birch-DEC) increased over time, except for the DEC in Block D in 1997. The decline in mean DBH values of the DEC trees was because several small trees attained a height of more than 1.3 m and were first measured in 1997 and a few of the larger DBH trees died. The mean DBH of living PL trees in Blocks C and D decreased from 2004 to 2013, due to the mortality of many of the larger PL trees caused by the mountain pine beetle (Dendroctonus ponderosae Hopkins) over this period. This mortality reduced the number of PL trees by 74% and 87% in Blocks C and D, respectively, by the end of 2013. The number of trees for the other species also declined over the 21-year measurement period since mortality considerably exceeded recruitment (ingrowth). The small amount of recruitment that occurred was generally only found in the thinned plots. Block D had a higher number of SX than Block C due to its proximity to the SBS zone; Block C had more DEC trees than Block D.

Choice of a Competition Index
Competition among trees reduces basal area growth and increases the probability of mortality. Researchers have used a number of spatially and non-spatially explicit indices to assess the effects of competition on individual tree growth [27][28][29]. We examined 13 widely used competition indices (five non-spatially explicit and eight spatially explicit indices) with respect to their ability to contribute to the prediction of tree basal area growth using a similar dataset to that employed in this study ( Table 2). Preliminary analyses using Equation (6) indicated that Martin and Ek's [30] competition index performed best overall for the various species groups, blocks, and growth periods based on fit statistics, and was consequently selected for use in this study (results not shown). Table 2. Sources of competition indices and corresponding formulas, tested for use in the individual tree basal area growth model. i, subject tree; j, neighbor tree; d i , DBH of subject tree (cm); d j , DBH of neighboring tree (cm); BA j , basal area of neighboring trees (cm 2 ); BAL i , basal area of trees larger than the subject tree (cm 2 plot −1 ); D q , quadratic mean diameter (cm); h i , height of subject tree (m); h j , height of neighboring tree (m); L i j , distance between subject tree and neighboring tree j; n, number of neighboring trees; mcw i , maximum crown width (m); a i j , area of influence zone overlap between neighboring tree and subject tree (m 2 ); A i , area of subject tree influence zone (m 2 ).

BA j
Lorimer [32] n j=1 d j d i Glover and Hool [33] d i D q Krajicek et al. [34] n i=1 π x mcw 2 i /4 Wykoff et al. [35]   Martin and Ek's [30] competition index is calculated as follows: where CI is the Martin and Ek competition index; d i is DBH of the subject tree (cm); d j is the DBH of the jth neighboring tree (cm); L ij is the distance between the subject tree and the neighboring tree; and n is the number of neighboring trees. The index relates the diameters d j of all potential competitors to the diameter of the subject tree d i , and then determines the sum of quotients to obtain the competition index. The quotients are weighted by the exponential function, where the contribution of a competitor to the competition index decreases with an increasing distance and decreasing diameter [42]. Neighbors were selected within a fixed 5 m radius of the subject tree, which a preliminary analysis indicated as the best search radius for identifying competitors based on their impact on the basal area growth of the subject trees. This search radius also matched with the buffer width surrounding each plot, allowing us to identify competitor trees for all the trees within each plot.

Basal Area Growth Models
Various functional forms have been used to model single-tree basal area increments (e.g., [43,44]). One commonly used nonlinear model is where PABAI is the periodic annual basal area increment (cm 2 yr −1 ) between two consecutive measurements; DBH is the diameter at breast height (cm) measured at the start of each growth period; β 1 , β 2 , and β 3 are the parameters to be estimated; and ε is the random error term. The term β 1 DBH β 2 indicates that basal area growth increases with the initial tree diameter, and exp (− ) is a term that is associated with a decline in the basal area increment as the size increases to prevent the unlimited growth of large-diameter trees [43,44]. Alternative functional forms for modeling basal area increments were also tested, including the models presented by [9,45]. Preliminary analyses of the alternative models showed that Equation (2) fitted our data well ( Figure A1) and, subsequently, it was chosen as the base model. Due to the complex structure and composition of tree species at the study site, a nonlinear mixed effects model [46] was fit for each main species using the package "nlme" in R 3.5.1 [47]. Nonlinear mixed models contain both fixed parameters common to all subjects, and random parameters specific to each subject. Following the single-level nonlinear mixed model function [46], Equation (2) can be expressed in general form as where PABAI ij is the periodic annual basal area increment, f is a nonlinear function with a plot-specific parameter vector x ij and the predictor vector v ij , ε ij is a normally distributed noise term, M is the number of plots, and n i is the number of observations on the ith plot. The plot-specific parameter vector x ij can be expressed as where β is a p-dimensional vector of fixed population parameters; b i is a q-dimensional random effects vector associated with the ith plot; A ij and B ij are design matrices for the fixed and random effects, respectively; and σ 2 D is a (general) variance covariance matrix. It is assumed that observations made on different plots are independent and that ε ij follow a normal distribution N 0, σ 2 and are independent of b i . When a new subject is available, the model can be calibrated for this subject by using information about the subject to estimate the empirical best linear unbiased predictors (EBLUPs) of the random parameters. A common option employed to predict random effects is the empirical Bayes approach [48], expressed asb whereb i is the estimated random effects parameter vector of EBLUPs,D is the estimated variancecovariance matrix of the random effects parameters, Z i is the estimated partial derivatives matrix with respect to the random effects parameters for the new observation,R i is the estimated variance-covariance matrix for the error term, andê i is the estimated residual vector. A random effect was initially specified for each of the parameters to account for the variability among plots. Preliminary fits based on the Akaike Information Criterion (AIC) values indicated that only the random effect associated with the β 2 parameter resulted in a smaller AIC. The random effects related to the β 1 and β 3 parameters were then dropped. The final model form was where b 1 is the random (plot) effect associated with parameter β 2 . The growth period was added as a covariate to Equation (6). Indicator (dummy) variables were used to represent the four growth periods and an initial analysis indicated that all parameters varied with growth period. The first growth period (1993-1997) was used as a reference. This model fitted our data well ( Figure A1) and, subsequently, it was chosen as the base model.
where β Gp is the parameter for the growth periods and GP is the indicator variable for the growth periods. Equation (7) was then extended to include an exponential function of the competition index to quantify the neighborhood effects on the growth of the subject tree: where β 4 is a parameter to be estimated, CI is Martin and Ek's competition index, and the other variables are as defined previously. Partitioning the competitive effects of neighbors into different species groups allowed us to incorporate potential differences in competitive relationships among species: where, β Fd , β Sx , β Pl , and β Dec are the parameters to be estimated for the competition index for Douglas-fir, spruce, lodgepole pine, and the broadleaf species, respectively (CI Fd , CI Sx , CI Pl , and CI Dec ), and the other variables are as defined previously. An initial analysis of the residuals showed an increase in variance as the DBH increased, violating the assumption of homogeneous variance. Modeling residual variance as a power function of the initial DBH best stabilized the residual variance, i.e., Throughout the model fitting process, the fit was examined using diagnostic plots of the standardized residuals versus the predicted values. The significance level was set at α = 0.05.

Analyses
We posited that the impact of competition on tree basal area growth could differ by species and block. We also wished to assess whether the impact of competition on basal area growth differed between small (defined here as trees from 2.0 to 9.9 cm DBH as of the 1993 measurements) and larger trees (defined here as trees ≥ 10.0 cm DBH as of the 1993 measurements) because of the possibility of competition effects varying with tree size [49]. We chose a 10 cm DBH as the cut point between the two size classes since many stand dynamics studies only begin to consider trees that have diameters larger than this value. We had sufficient numbers of trees in the dataset to fit stable PABAI models (i.e., Equations (7)-(9) for the following groups of subject trees: (1) small FD in Block C; (2) large FD in Block C; (3) small FD in Block D; (4) large FD in Block D; (5) large SX in Block D; and (6) large PL in Block D. Due to mortality caused by a mountain pine beetle infestation, we only had sufficient PL to fit PABAI models for the first two growth periods (i.e., up to 2004). We were able to fit PABAI models for all four growth periods for all other groups of subject trees. As subject trees, we only used trees that were present in 1993 and survived until 2013. However, all trees that were alive during a particular measurement interval were considered potential competitors with a subject tree during that measurement interval. The DEC trees were considered competitors, but were not modeled as subject trees.
We used the likelihood ratio test (LRT) to determine the significance of the improvement of the extended model (Equation (8)) over the base model (Equation (7)) for each group of subject trees. We then compared, for each group, the extended model (Equation (8)) and the partitioned CI model (Equation (9)) using AIC to determine whether identifying the species of the competitor allowed better predictions of PABAI. The goodness-of-fit for each model was assessed with the root mean square error (RMSE); the marginal R 2 , which is a measure of the proportion of variance explained by the fixed factors alone; and conditional R 2 , which is the variance explained by both fixed and random factors. The significance of the indicator variable (growth periods) was assessed using an F test.
We assessed the growth of different-sized trees for FD, SX, and PL under no competition and then under different levels of generic competition (the 10th, 50th, and 90th percentiles of the Martin and Ek competition index for a particular measurement period) using Equation (8). For these model applications, we only used the fixed-effect parameter estimates (i.e., population-averaged predictions) [46,49].

Change in Stand Characteristics
The mean DBH for small FD changed from 6.4 to 8.3 cm in Block C and 6.4 to 8.1 cm in Block D, while the mean DBH for larger FD changed from 14.9 to 20.1 cm in Block C and 14.4 to 19.0 cm in Block D over the growth periods ( Table 3). The mean DBH of large SX trees changed from 15.6 to 18.5 cm over the growth periods, while that of PL increased from 16.8 to 19.1 cm over the first two growth periods. PABAI for all small and large trees, irrespective of species, decreased in the second growth period (1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004), but increased in the third growth period (2004)(2005)(2006)(2007)(2008), and declined again in the last growth period in Block D. The pattern in Block C was similar, except that there was a slight decline in PABAI for small FD in the third growth period. Table 3. Means and standard deviations (in brackets) for surviving trees from 1993 to 2013 by tree size and species. DBH is the diameter at breast height (cm); BA is the tree basal area (cm 2 ); PABAI is the periodic annual basal area increment (cm 2 yr −1 ); n is the number of trees; FD is interior Douglas-fir; SX is interior spruce; and PL is lodgepole pine. Small trees were between 2.0 and 9.9 cm at the time of the 1993 measurement and large trees were ≥ 10.0 cm DBH at the 1993 measurement.

Models for Predicting PABAI
The likelihood ratio test (LRT) indicated that using a single CI (Equation (8), Figure A3) significantly improved the fit over the base model (Equation (7), Figure A2) for all species and blocks ( Table 4). All parameter estimates (coefficients) associated with the Martin and Ek (30) competition index were negative and all were significant for the subject trees in Blocks C and D (Figure 2a,b) for all growth periods. The impact of competition on PABAI varied for small and large FD across the growth periods. For Block C, there was a significant effect of competition between growth periods for small FD (F 1,3 = 8.08, p < 0.0001) and large FD (F 1,3 = 8.28, p < 0.0001). Competition was more pronounced in the first growth period and less in the last growth period for both small and large FD. For Block D, the competition index coefficients were similar for small and large FD in growth periods 1, 2, and 4; however, in the third growth period, the magnitude of the coefficient was smaller for the small FD. Large SX in Block D had the highest magnitude competition index coefficient of all the conifers in the third growth period, but it was the lowest in the second growth period.  (9) (the model with a partitioned competition index). FD is interior Douglas-fir; SX is interior spruce; PL is lodgepole pine; "small" are trees from 2.0 to 9.9 cm DBH in 1993; "large" are trees ≥ 10.0 cm DBH in 1993; LRT is the statistic for the likelihood ratio test for comparing the base model and the single CI model; AIC is the Akaike Information Criterion used for comparing the single CI model to the partitioned CI model; RMSE is the root mean square error; Marginal and Conditional R 2 are the pseudo R 2 for the fixed factor alone, and both fixed and random factors, respectively; and * indicates where the LRT was significant at α = 0.05. The smallest AIC and RMSE for a species group and block are shown in bold.

Blocks
Species Tree Size Equation (  Partitioning the competition index by species (Equation (9), Figure 3) provided a better fit than the model that used only a single competition index (Equation (8), Figure A3) for small FD in Block C and large FD in Block D (Table 4). We did not find any gain when predicting PABAI using Equation (9) for large FD in Block C, and small FD, large SX, and PL in Block D. All models showed a pseudo R 2 between 53% and 80% when both fixed and random factors were included and between 37% and 75% when only fixed factors were considered. The worst fit was obtained for lodgepole pine, where the pseudo R 2 ranged between 53% and 61% and 37% and 54%, including all or only fixed factors, respectively.  Partitioning the competition index by species (Equation 9, Figure 3) provided a better fit than the model that used only a single competition index (Equation 8, Appendix A3) for small FD in Block C and large FD in Block D (Table 4). We did not find any gain when predicting PABAI using Equation 9 for large FD in Block C, and small FD, large SX, and PL in Block D. All models showed a pseudo between 53% and 80% when both fixed and random factors were included and between 37% and 75% when only fixed factors were considered. The worst fit was obtained for lodgepole pine, where the pseudo ranged between 53% and 61% and 37% and 54%, including all or only fixed factors, respectively. We interpreted the competitive relationships among species by looking at the relative magnitude of the coefficients for the partitioned competition index in Equation 9 for the four growth periods. The coefficient for each subsequent growth period (second, third, and fourth) was added to the coefficient of the reference coefficient (for the first growth period) to obtain the magnitude coefficient for the specified period. For small FD in Block C, the parameter estimates for FD competition were more negative than those of SX, PL, and DEC competition for all but the 2008-2013 growth period, indicating that competition from neighboring FD (i.e., intraspecific competition) was heavier than competition from the other species (Table 5). This pattern was also similar for large FD in Block C, where parameter estimates for FD competition were more negative than those of SX, PL, and DEC competition in the first and third growth periods. There was a significant effect of FD competition between growth periods for small and large FD in Block C, but not small and large FD in Block D (Table 6). For both large and small FD subject trees in Block D, the parameter estimates for We interpreted the competitive relationships among species by looking at the relative magnitude of the coefficients for the partitioned competition index in Equation (9) for the four growth periods. The coefficient for each subsequent growth period (second, third, and fourth) was added to the coefficient of the reference coefficient (for the first growth period) to obtain the magnitude coefficient for the specified period. For small FD in Block C, the parameter estimates for FD competition were more negative than those of SX, PL, and DEC competition for all but the 2008-2013 growth period, indicating that competition from neighboring FD (i.e., intraspecific competition) was heavier than competition from the other species (Table 5). This pattern was also similar for large FD in Block C, where parameter estimates for FD competition were more negative than those of SX, PL, and DEC competition in the first and third growth periods. There was a significant effect of FD competition between growth periods for small and large FD in Block C, but not small and large FD in Block D (Table 6). For both large and small FD subject trees in Block D, the parameter estimates for FD competition were not the most negative, indicating that interspecific competition was higher than intraspecific competition, except for large FD in the third growth period (Table 5). PL and SX were the heaviest competitors for large FD in two and one of the four growth periods, respectively. For small FD, PL was the heaviest competitor for two of the four growth periods, and SX and DEC were the heaviest competitors in the second and third growth periods, respectively. For large SX subject trees, the heaviest competition was from PL in three of the four growth periods and from DEC in the other period. For large PL subject trees for the two growth periods for which models were fit, the heaviest competition was also from DEC in the first growth period and neighboring PL in the second growth period. Table 5. Estimated parameters of the periodic annual basal area growth models (Equation (9) Table 6. F tests for the significance of the fixed effects estimated for the partitioned competition index (Equation (9)) for Blocks C and D. The fixed effect parameter estimates are presented in Table 5

Effect of Competition on the Growth of the Subject Trees
In general, FD was projected to have higher PABAI than PL and SX for all growth periods and, not surprisingly, the growth rate declined with increasing competition levels for the three species (Figure 4a-d). Comparing the potential BA (i.e., the predicted basal area when a tree is free from competition- Figure 4a) to when a tree of a similar size is experiencing high competition (90th percentile of the competition index- Figure 4d), the loss in BA increased over time. For example, for a tree of a 20 cm DBH in 1993 in Block D, the decline in potential basal area with increasing levels of competition ranged from 10% to 39% for FD and 6% to 34% for SX over the four growth periods. The BA for PL declined by 11-15% over the first and second growth periods. Large FD were predicted to grow slightly better in Block C compared to Block D at lower levels of competition; however, the growth of FD in Block D declined less than that in Block C at higher levels of competition, resulting in similar or slightly better growth in Block D (Figure 5a-d). Similarly, small FD of a 5 or 8 cm DBH in 1993 grew slightly faster in Block C than Block D at low levels of competition, but the predicted growth was quite similar between the blocks for higher levels of competition (Figure 6a,d). FD that were 2 cm in 1993 were predicted to grow similarly in both blocks at all levels of competition. growth was quite similar between the blocks for higher levels of competition (Figures 6a and 6d). FD that were 2 cm in 1993 were predicted to grow similarly in both blocks at all levels of competition.
(a) (b) (c) (d)    growth was quite similar between the blocks for higher levels of competition (Figures 6a and 6d). FD that were 2 cm in 1993 were predicted to grow similarly in both blocks at all levels of competition.

Discussion
There is a finite amount of total growing space (considered as the sum total of light, moisture, and nutrients) available on a site to support trees. If there are a sufficient number of trees to occupy that growing space, then trees compete for those resources [4]. If competitive reduction exists among

Discussion
There is a finite amount of total growing space (considered as the sum total of light, moisture, and nutrients) available on a site to support trees. If there are a sufficient number of trees to occupy that growing space, then trees compete for those resources [4]. If competitive reduction exists among tree species in a mixed-species stand, it suggests that there would be differences among species in terms of their competitive impact on each other's growth. In other words, one might expect the growth of a tree of a given species to be impacted more intensively by trees of the same species than those of other species, all other factors being equal. On the other hand, if a species' growth is impacted more strongly by competition from trees of other species than of the same species, it would suggest that a process of competitive exclusion is taking place. It is quite possible that both processes are taking place within naturally established mixed-species stands, such as the stands we examined in this study.
We investigated whether neighboring trees compete with a subject tree independently of species identity and examined the differences in basal area growth among FD, SX, and PL over time under different degrees of competition. The Martin and Ek [30] competition index, found to be suitable for this study, is characterized by an asymmetrical relationship in which competitors share limited resources disproportionately in relation to their relative sizes [50]. This index is thought to reflect aboveground competition for light in which larger competitors pre-empt resources, making them unavailable to smaller individuals. Although multiple types of interactions can occur simultaneously in mixtures, our study found asymmetric competition to be more relevant than symmetric competition, suggesting that the main driver limiting growth in these stands is aboveground competition for light. This is not surprising as the DEC species and PL are intolerant of shade, and considered pioneer species [51], and SX and FD are considered to be moderately shade-tolerant [51,52]. This agrees with previous studies that highlight the importance of light when comparing size-asymmetric and size-symmetric competition in mixed-species stands [45,53,54].
Interactions with neighbors exerted a negative influence on basal area growth for all the subject trees we considered and these influences varied across time. These fluctuations in the intensity of competition on basal area growth over time show that temporal dynamics in resource availability, likely caused by year-to-year variations in weather (particularly precipitation patterns) and mortality of many of the larger PL, influenced the magnitude of interactions between species. Although other interactions such as facilitation and competitive reduction among species have been observed in mixed-species forests by other studies [15,55], this was not found in our study, except for small FD in Block C and large FD in Block D. In these instances, the model with the competition index partitioned by species (Equation (9)) fit best. Block C is slightly drier than Block D, and it is possible that the smaller FD is better able to compete for moisture with the other species groups than it is able to do so with itself. We did not see this phenomenon for small FD on Block D, and in fact, the model with a single competition index (Equation (8)) generally fitted best for the subject trees on that block. Further, the competition indices coefficients for the various Equation (9) fits in Block D indicated that species other than FD were exerting the most influence on the basal area growth of all the subject tree species groups. This might indicate that in the moister conditions found at Block D, species other than FD are better able to exploit belowground resources. This is consistent with this block location being considered transitional from the IDF BEC zone to the SBS BEC zone, where SX is considered to be the climax species.
The mixed results we found with respect to the impact of species identity on predicting basal area growth are consistent with what has been shown across other studies. Some studies have demonstrated differences in the species identity of neighbors in some species and the functional equivalence of neighbors in others (e.g., [56][57][58]). The apparent equivalence of neighbors in mixed-species stands has been attributed to a small sample size, which made it difficult for models to detect interspecific differences in competitive effects [20,57,58]. However, in our study, we had sufficient numbers of trees both as subject trees and as competitors. Instead, our results may be due to the complex processes occurring in multi-aged mixed-species forests where relationships between competition and tree growth may be diverse, with different species and tree sizes varying in their competitive interactions with neighboring species. In addition, the complexity of the competition environment at our sites could be compounded by variability in stochastic disturbance events or partial harvesting, which resulted in clusters of regeneration following disturbance. For example, LeMay et al. [59], in similar, but non-harvested, stands to those we examined in this study, and Druckenbrod et al. [60], found the clustering of trees at shorter distances.
The differentiation of intraspecific and interspecific competition provides useful information for understanding factors influencing tree growth and interactions between species [61]. In our study, for small FD in Block C, where a better fit was obtained when the competition index was partitioned into species or species groups, there were differences in species-specific effects. FD had, on average, stronger competitive impacts from intraspecific neighbors than from other species. Similarly, large FD in Block C and PL in Block D had stronger competition from intraspecific neighbors. Many other studies of mixed-species stands have found that intraspecific neighbors are stronger competitors than interspecific neighbors [15,45,62]. However, this effect can be positive or negative, depending on the species characteristics and site conditions, such as shade tolerance [63], stand density, and stand development [64].
We observed a higher PABAI for FD than for PL and SX, irrespective of the competition intensity. Similar growth patterns have been reported on naturally regenerated saplings of PL and FD along a wide range of light conditions, where PL showed a greater reduction in lateral growth than FD with a decreasing light availability [51]. Eis et al. [65] found that SX grows slowly, even if it started growing at the same time as PL. We expected small and large FD in Block C to grow better than small and large FD in Block D; however, at a high competition intensity, FD trees in Block D were growing at the same rate or slightly better than FD trees in Block C. This may be due to the nature of the interspecific competition experienced by FD trees in Block D versus the primarily intraspecific competition in Block C. Also, Block D is moister than Block C and belowground competition for moisture may be more limiting to growth at higher competition levels in Block C.

Conclusions
There was little differentiation between intraspecific and interspecific competition for the large subject trees we observed, which is consistent with competition for light in situations where the two primary species (FD and SX) are similar in shade tolerance (mid-tolerant). PL is shade-intolerant, but it grew either in less dense pockets in the stands or was one of the taller trees. Again, the species of its competitors appeared to have little impact on its basal area growth. In contrast, smaller FD were affected more heavily by intraspecific competition, suggesting that belowground competition with other FD, likely for moisture, could be playing a role, along with competition for light, in their growth. We did not have sufficient numbers of small SX trees to fit stable basal area growth prediction equations for that stand component, so we were not able to determine whether the basal area growth of these trees was also more heavily affected by competition with FD. Our models suggested that FD basal area growth, for a given tree size (DBH), was higher than that of PL and SX, irrespective of the competition level and block.
It is clear that competition is a complex process in mixed-species mixed-sized natural stands and is undoubtedly affected by a number of factors, including tree sizes, species composition, and developmental history. The relationships between periodic basal area growth and competition varied with species identity over time, likely as a result of both variations in weather (that affects moisture availability during the growing season) and disturbances. Our results provide an understanding of how IDF stands will develop over time and information on species interactions that could help forest managers design more effective silvicultural prescriptions.