Abundance and Impacts of Competing Species on Conifer Regeneration Following Careful Logging in the Eastern Canadian Boreal Forest

Managing competing vegetation is crucial in stand establishment strategies; forecasting the abundance, composition, and impact of competing vegetation after harvesting is needed to optimize silviculture scenarios and maintain long-term site productivity. Our main objective was to identify factors influencing the short-term abundance and composition of competing vegetation over a large area of the Canadian boreal forest. Our second objective was to better understand the mid-term evolution of the regeneration/competing vegetation complex in cases of marginal regeneration conditions. We used operational regeneration surveys of 4471 transects sampled ≈5 years after harvesting that contained data on regeneration, competing vegetation, elevation, ecological classification, soil attributes, and pre-harvest forest stands. We performed a redundancy analysis to identify the relationships between competing vegetation, harvesting and biophysical variables. We then estimated the probability of observing a given competing species cover based on these variables. In 2015, we re-sampled a portion of the sites, where conifer regeneration was marginal early after harvesting, to assess the temporal impact of different competing levels and species groups on the free-to-grow stocking, vigour and basal area of softwood regeneration. Results from the first inventory showed that, after careful logging around advance growth, ericaceous shrubs and hardwoods were not associated with the same sets of site attributes. Ericaceous shrubs were mainly found on low fertility sites associated with black spruce (Picea mariana (Mill.) BSP) or jack pine (Pinus banksiana Lamb.). The distinction between suitable environments for commercial shade-intolerant hardwoods and non-commercial hardwoods was less clear, as they responded similarly to many variables. Analysis of data from the second inventory showed a significant improvement in conifer free-to-grow stocking when commercial shade-intolerant hardwood competing levels were low (stocking 0%–40%) and when ericaceous shrubs competing levels were moderate (percent cover 26%–75%). In these conditions of marginal regeneration, the different types and intensities of competition did not affect the vigour or basal area of softwood regeneration, 9–14 years after harvesting.


Introduction
Stand renewal is a crucial step in any silvicultural system.At this stage, desired species can experience competition by a suite of species that can impact longer-term stand composition and though they have commercial values, these hardwoods are often considered as competing species because they can impede the growth and survival of the more economically valuable conifers.Other common competing hardwoods are mountain maple (Acer spicatum Lam.), pin cherry (Prunus pensylvanica L. f.) and speckled alder (Alnus incana subsp.rugosa (Du Roi) J. Clausen) [18].Ericaceous competition is mainly caused by bog Labrador tea (Rhododendron groenlandicum (Oeder) K.A. Kron and Judd), sheep laurel (Kalmia angustifolia L.) and blueberries (Vaccinium spp.) [19].

Post-Harvest Inventory
To understand factors influencing the short-term abundance and composition of competing vegetation, we used data collected operationally between 2007 and 2011 by AbitibiBowater Inc. Ecological regions are those defined by [16].

Post-Harvest Inventory
To understand factors influencing the short-term abundance and composition of competing vegetation, we used data collected operationally between 2007 and 2011 by AbitibiBowater Inc. (now Resolute Forest Products Inc.) in post-harvest regenerating stands.Monitoring was conducted within the first eight years (usually four or five years) following operational careful logging around advance growth.Sampling was conducted in north-oriented transects of 10 micro-plots separated 5 m from each other (Figure 2a).Each micro-plot consisted in two superimposed circular plots with radii of either 1.13 m (plots ≈ 4 m 2 ) for conifer observations, or 1.69 m (plots ≈ 9 m 2 ) for hardwood observations [20].Data collected in each micro-plot consisted in (1) presence/absence of regenerating trees by species and height class (15-30 cm, 30-60 cm, 60-100 cm, >100 cm), and (2) percent cover for groups of competing species (ericaceous shrubs and non-commercial hardwoods) using 25% classes (0%, 1%-25%, 26%-50%, 51%-75%, 76%-100%).Percent cover was collected by groups for non-commercial hardwoods and ericaceous shrubs.Using data from the micro-plots, we calculated stocking at the transect level by species (spruce sp., balsam fir, jack pine, trembling aspen, paper birch) and species groups (conifers, hardwoods).
conducted within the first eight years (usually four or five years) following operational careful logging around advance growth.Sampling was conducted in north-oriented transects of 10 microplots separated 5 m from each other (Figure 2a).Each micro-plot consisted in two superimposed circular plots with radii of either 1.13 m (plots ≈ 4 m 2 ) for conifer observations, or 1.69 m (plots ≈ 9 m 2 ) for hardwood observations [20].Data collected in each micro-plot consisted in (1) presence/absence of regenerating trees by species and height class (15-30 cm, 30-60 cm, 60-100 cm, >100 cm), and (2) percent cover for groups of competing species (ericaceous shrubs and noncommercial hardwoods) using 25% classes (0%, 1%-25%, 26%-50%, 51%-75%, 76%-100%).Percent cover was collected by groups for non-commercial hardwoods and ericaceous shrubs.Using data from the micro-plots, we calculated stocking at the transect level by species (spruce sp., balsam fir, jack pine, trembling aspen, paper birch) and species groups (conifers, hardwoods).For each transect, metadata from the inventory was cross-validated using external data sources.Time since harvesting was determined from known harvesting dates, and harvesting seasons were verified using Landsat images.The presence/absence of delimbing areas, skidding cones and orthogonal trails were used to validate the use of cut-to-length or full tree harvesting methods.Sites that had been affected by wildfires since 1969 were eliminated from the database.
Transect altitude was determined from Quebec's topographic database.Pre-harvest forest stands and soil attributes were derived from the governmental forest ecological survey dataset [21].Sites classified as alder groves, dry barrens or wet barrens were discarded, being considered as unproductive from a forest management perspective.As a result, 4471 transects were retained, representing 43,972 microplots distributed in 1277 cutblocks.

Re-Sampling of Sites with Marginal Regeneration Conditions
During summer 2015, we re-sampled 72 of the transects from the previous inventory to evaluate the evolution of regeneration in stands with free-to-grow stockings between 40% and 60% that were still accessible.Seedlings were considered free-to-grow when no vegetation within a 1 m radius exceeded half the height of the target tree [22].The choice of focusing on this range of stocking was based on the fact that, five years after harvesting, a conifer stocking of 60% is about the minimum threshold for a stand to reach maximum yield at maturity [23] and a site is generally considered as non-regenerated when its commercial species stocking is less than 40% [24].Between For each transect, metadata from the inventory was cross-validated using external data sources.Time since harvesting was determined from known harvesting dates, and harvesting seasons were verified using Landsat images.The presence/absence of delimbing areas, skidding cones and orthogonal trails were used to validate the use of cut-to-length or full tree harvesting methods.Sites that had been affected by wildfires since 1969 were eliminated from the database.
Transect altitude was determined from Quebec's topographic database.Pre-harvest forest stands and soil attributes were derived from the governmental forest ecological survey dataset [21].Sites classified as alder groves, dry barrens or wet barrens were discarded, being considered as unproductive from a forest management perspective.As a result, 4471 transects were retained, representing 43,972 microplots distributed in 1277 cutblocks.

Re-Sampling of Sites with Marginal Regeneration Conditions
During summer 2015, we re-sampled 72 of the transects from the previous inventory to evaluate the evolution of regeneration in stands with free-to-grow stockings between 40% and 60% that were still accessible.Seedlings were considered free-to-grow when no vegetation within a 1 m radius exceeded half the height of the target tree [22].The choice of focusing on this range of stocking was based on the fact that, five years after harvesting, a conifer stocking of 60% is about the minimum threshold for a stand to reach maximum yield at maturity [23] and a site is generally considered as non-regenerated when its commercial species stocking is less than 40% [24].Between these two thresholds, conifer establishment success is uncertain.Sites were divided into three strata representing the initial level of competing vegetation (low, medium, high; see Table 1) and two strata representing the initial type of dominant competing vegetation (Table 1).Sites initially dominated by non-commercial hardwoods were not included in the stratification due to a lack of transects belonging to the "high" level of competition.Transects that had been submitted to mechanical site preparation and/or plantation were also avoided.
Based on the initial sampling, we established transects of 10 circular micro-plots separated by 5 m from each other, except for the fifth whose centre was distanced 7 m from the adjacent micro-plots (Figure 2b).Micro-plot areas and vegetation cover classes were similar to the initial inventory.In each micro-plot, we noted the presence/absence and free-to-grow status of regeneration by species and height class.Saplings were counted in a 5.64 m radius plot superimposed over the fifth micro-plot.In this micro-plot, three sample trees representative of the modal height of the regeneration were selected; we noted the species, height, length of the live crown, length of the terminal shoot and length of the longest lateral branch of the last whorl.Tree height and length of the live crown were used to calculate live crown ratio [25].Length of the terminal shoot and length of the longest lateral branch of the last whorl were used to calculate apical dominance ratio [26].Live crown ratio and apical dominance ratio are two indices used to estimate conifer vigour.Ericaceous shrubs: low = percent cover between 1%-25%, medium = percent cover between 26%-75%, high = percent cover between 76%-100%.Commercial shade-intolerant hardwoods: low = stocking between 0%-40%, medium = stocking between 40%-60%, high = stocking between 60%-100%.

Data Analysis
To identify factors influencing the abundance and composition of competing vegetation, we performed a redundancy analysis [27] (RDA) on data from the first inventory, using the vegan package [28] in R version 3.3.3[29].We split data into two matrices.The first matrix (the Y matrix) contained the dependent variables, namely, the percent cover of ericaceous shrubs, the percent cover of non-commercial hardwoods and the stocking of commercial shade-intolerant hardwoods.The second matrix (the X matrix) contained explanatory variables, namely, harvesting and forest characteristics such as pre-harvest stand dominant species group, type and year of the original disturbances, drainage, surficial materials, ecological region, elevation, harvesting method and harvesting season.Before the RDA, we kept explanatory variables containing more than 20 observations (n > 20), to detect reasonable size effects with reasonable power [30].We also applied a Hellinger transformation to the Y matrix to give low weights to variables with low counts and many zeros and thus maintain an ecologically meaningful distance among sites in the ordination [31].We also performed a forward selection using a double-stopping criterion [32] to select the best subset of explanatory variables and avoid strong multicollinearity.We verified multicollinearity by making sure variance inflation factors (VIF) were <10 for the explanatory variables [27].
To estimate the probability of observing a given competing species group according to harvesting and environmental characteristics, we performed a linear mixed-effects analysis [33] (LMM) of the commercial shade-intolerant hardwood stocking with the lme4 package [34] in R. Degrees of freedom and p-values for the LMM were obtained using the lmerTest package [35].We also performed two cumulative link mixed analyses [36] (CLMM) of the cover of non-commercial hardwoods and ericaceous shrubs with the ordinal package [37].For LMM and CLMMs, cutblocks were used as a random effect factor and the X variables identified in the RDA were used as fixed effects.Commercial shade-intolerant hardwood stocking was used as a dependent variable for the LMM and percent cover of ericaceous shrubs or commercial shade-intolerant hardwoods were used as dependent variables for the CLMMs.We calculated a pseudo-R 2 for the LMM with the MuMIn package [38].
To evaluate the evolution of regeneration in stands with free-to-grow stockings between 40% and 60%, we performed four mixed analyses of variance (ANOVAs) [39] (also referred to as "nested ANOVA") using the nlme package [40] and the anova function of R. When a significant effect was detected (p < 0.05) for an interaction or a main effect, we proceeded to a pairwise comparison of the least square means with the lsmeans package [41].For the first mixed ANOVA, we checked whether the competing cover type and level affected free-to-grow conifer stocking.To do so, we used the Forests 2019, 10, 177 6 of 17 conifer stocking measured in both inventories as dependent variables and the six strata of competing cover (Table 1), as well as the time of the inventory (1 for the first inventory, 2 for the second) as fixed effects.For the random effects, we used sites (to account for the repeated measures) and cutblocks (to account for the nested design).We tested for interaction between time and competing vegetation cover.For the second and third analyses, we asked if the competing level measured during the first inventory influenced the vigour of the conifer regeneration assessed at the second inventory.To do so, we relied on the two vigour indices (live crown ratio and apical index) and performed a mixed ANOVA using each of these as dependent variables and strata of competing covers as a fixed effect.We used cutblocks as a random effect to account for the nested design.Finally, we asked if the competing level measured in the first inventory influenced sapling basal areas measured in the second inventory.To do so, we used sapling basal area as a dependent variable and the six strata of competing covers as a fixed effect.To account for the nested design, we used cutblocks as a random effect.

Variables Significantly Linked to Competing Vegetation Cover
Explanatory variables (Table 2) selected for the RDA accounted for 35.3% of the variation in the composition of competing vegetation; axes 1 and 2, respectively, explained 32.98% and 2.02% of the variance (Figure 3).Permutation tests indicated that the global model and canonical axes were significant at p = 0.001.Axis 1, which contains most of the information, contrasts conditions associated with ericaceous shrubs and those associated with other competing vegetation types.The ericaceous shrubs were positively correlated with imperfect drainage, altitude, hills of Lake Péribonka and pre-harvest stands dominated by black spruce (Figure 3).Hardwoods were negatively correlated with all the above and positively correlated to the hills surrounding Lake Saint-Jean, tills 50-100 cm thick, pre-harvest stands dominated by paper birch and pre-harvest stands dominated by both balsam fir and paper birch.Non-commercial hardwoods were also positively correlated to partial harvesting.competing species groups (black arrows) and explanatory variables selected using a forward selection approach (grey arrows).Only variables with the highest contribution to axes RDA1 or RDA2 are labeled (i.e., coordinate on one axis was >90th quantile or <10th quantile of the distribution of variables coordinates on the same axis).Refer to Table 2 for variables description.
The linear mixed-effects (LMM) analysis predicting the probability of observing a given level of commercial shade-intolerant hardwoods had a marginal pseudo-R 2 of 0.34.Its coefficient estimates (Table 3) show that the chances of observing this group significantly increased in the presence of good drainage, full tree-harvesting, tills 50-100 cm thick, and pre-harvest stands dominated by paper birch and accompanied by trembling aspen.These chances significantly decreased in the presence of imperfect drainage, ecological regions 6g or 6h, winter harvesting, pre-harvest stands originating from total windthrow, clearcutting or mild epidemic, pre-harvest paper birch stands, pre-harvest stands containing black spruce or jack pine as a dominant or codominant species, disintegration moraine, ice-contact deposit, and outwash.species groups (black arrows) and explanatory variables selected using a forward selection approach (grey arrows).Only variables with the highest contribution to axes RDA1 or RDA2 are labeled (i.e., coordinate on one axis was >90th quantile or <10th quantile of the distribution of variables coordinates on the same axis).Refer to Table 2 for variables description.Thin organic deposits a of the hardwood basal area; b of the conifer basal area.
The linear mixed-effects (LMM) analysis predicting the probability of observing a given level of commercial shade-intolerant hardwoods had a marginal pseudo-R 2 of 0.34.Its coefficient estimates (Table 3) show that the chances of observing this group significantly increased in the presence of good drainage, full tree-harvesting, tills 50-100 cm thick, and pre-harvest stands dominated by paper birch and accompanied by trembling aspen.These chances significantly decreased in the presence of imperfect drainage, ecological regions 6g or 6h, winter harvesting, pre-harvest stands originating from total windthrow, clearcutting or mild epidemic, pre-harvest paper birch stands, pre-harvest stands containing black spruce or jack pine as a dominant or co-dominant species, disintegration moraine, ice-contact deposit, and outwash.Proportional odds assumption was not met for the CLMMs of the ericaceous shrubs and non-commercial hardwoods; the explanatory variables did not have the same effect on the odds from one threshold to another.This assumption is, however, rarely met [30], especially in the presence of many explanatory variables and a large sample size [42], as in our case.However, these models can still be useful; estimates obtained from both models provided a general idea of the changes in competing species' percent cover induced by the explanatory variables.Furthermore, for both CLMMs, altitude was excluded to avoid convergence problems.
The CLMM's coefficient estimates (Table 4) show that the chances of observing non-commercial hardwoods significantly increased in the presence of good drainage, in ecological regions 4e or 5d, in pre-harvest stands originating from partial windthrow or partial harvesting, in pre-harvest paper birch stands and pre-harvest stands dominated by paper birch with balsam fir or trembling aspen Forests 2019, 10, 177 9 of 17 as companion species.These chances significantly decreased in the presence of imperfect or poor drainage, in ecological region 6h, in pre-harvest stands originating from total windthrow, in pre-harvest stands dominated by black spruce, in pre-harvest stands dominated by jack pine, in pre-harvest stands dominated by balsam fir and accompanied by black spruce, in disintegration moraine, ice-contact deposit and outwash.The chances of observing ericaceous shrubs (Table 5) significantly increased in ecological region 6h, after winter harvesting, in pre-harvest stands originating from total windthrow, clearcutting, Forests 2019, 10, 177 10 of 17 in pre-harvest stands dominated by paper birch and accompanied by black spruce, pre-harvest stands dominated by black spruce or jack pine, in ice-contact deposits and in outwash.Those chances significantly decreased in the ecological regions 4e or 5d where there was full tree harvesting, cut-to-length logging, pre-harvest stands originating from partial windthrow or mild epidemic, pre-harvest paper birch stands, pre-harvest stands dominated by paper birch with balsam fir or trembling aspen as companion species, pre-harvest stands dominated by balsam fir with paper birch as companion species.

Temporal Impact of Competing Vegetation in Conditions of Marginal Regeneration
There was a significant interaction between time and competing cover strata (F(5, 66) = 2.7, p = 0.028).However, pairwise comparisons of the least square means showed that the increase in stocking was only statistically significant for low competing levels of commercial shade-intolerant hardwoods (t(66) = −5.992,p < 0.001) and moderate competing levels of ericaceous shrubs (t(66) = −5.653,p < 0.001) (Figure 4).We did not detect statistically significant differences over time for low competing levels of ericaceous shrubs (t(66) = −3.205,p < 0.081), high competing levels of ericaceous shrubs (t(66) = −2.126,p < 0.607), moderate competing levels of commercial shade-intolerant hardwoods (t(66) = −3.300,p < 0.063) and high competing levels of commercial shade-intolerant hardwoods (t(66) = −2.709,p < 0.245) (Figure 4).Free-to-grow conifer stocking had a tendency to improve over time in every strata of competing cover (Figure 4).On average, at the time of the second inventory, we observed free-to-grow conifer stocking exceeding 60% for every strata (Figure 4).Refer to Table 2 for variable description.Bold indicate significance at α = 0.05.

Temporal Impact of Competing Vegetation in Conditions of Marginal Regeneration
There was a significant interaction between time and competing cover strata (F(5, 66) = 2.7, p = 0.028).However, pairwise comparisons of the least square means showed that the increase in stocking was only statistically significant for low competing levels of commercial shade-intolerant hardwoods (t(66) = -5.992,p < 0.001) and moderate competing levels of ericaceous shrubs (t(66) = -5.653,p < 0.001) (Figure 4).We did not detect statistically significant differences over time for low competing levels of ericaceous shrubs (t(66) = -3.205,p < 0.081), high competing levels of ericaceous shrubs (t(66) = -2.126,p < 0.607), moderate competing levels of commercial shade-intolerant hardwoods (t(66) = -3.300,p < 0.063) and high competing levels of commercial shade-intolerant hardwoods (t(66) = -2.709,p < 0.245) (Figure 4).Free-to-grow conifer stocking had a tendency to improve over time in every strata of competing cover (Figure 4).On average, at the time of the second inventory, we observed free-to-grow conifer stocking exceeding 60% for every strata (Figure 4).There were no significant difference in live crown ratio (F(5, 13) = 2.296, p = 0.106) and apical index (F(5, 13) = 1.030, p = 0.440) according to competing cover strata.Most sampled trees were vigorous.For each competing cover strata, average live crown ratio was higher than 70% and average apical index was higher than 1.20.Least square means did not detect any significant difference in sapling basal area between the different competing cover strata.Average sapling basal area ranged from 1.67 to 5.19 m 2 /ha for each competing cover strata.

Variables Significantly Linked to Competing Vegetation Cover
Our results contrast the nature of competing vegetation between two major ecosystem types differing in terms of climate and fertility.Ericaceous shrubs were more commonly associated with black spruce or jack pine ecosystems while the opposite was true for intolerant hardwoods and non-commercial hardwoods.There were no significant difference in live crown ratio (F(5, 13) = 2.296, p = 0.106) and apical index (F(5, 13) = 1.030, p = 0.440) according to competing cover strata.Most sampled trees were vigorous.For each competing cover strata, average live crown ratio was higher than 70% and average apical index was higher than 1.20.Least square means did not detect any significant difference in sapling basal area between the different competing cover strata.Average sapling basal area ranged from 1.67 to 5.19 m 2 /ha for each competing cover strata.

Variables Significantly Linked to Competing Vegetation Cover
Our results contrast the nature of competing vegetation between two major ecosystem types differing in terms of climate and fertility.Ericaceous shrubs were more commonly associated with black spruce or jack pine ecosystems while the opposite was true for intolerant hardwoods and non-commercial hardwoods.
Ericaceous shrubs were more associated with colder ecological regions located to the north of the study area or sites located at higher altitudes.They were also associated with black spruce and jack pine stands that are typical in these conditions [43,44].In these stands, ericaceous shrubs are often present in the understory and can expand after harvesting [14].Cold climate slows down decomposition rates and favours the accumulation of organic matter [45], a substrate more favourable to black spruce and ericaceous shrub growth.These species also often occur together on poorly drained soils [19].Poor drainage favours the accumulation of organic matter, hence reducing nutrient availability, a condition less favourable to hardwood development [46].The abundance of sphagnum mosses on poorly drained soils favours the creation of acidic, wet and cold soils that decrease the decomposition rate of the organic matter, nutrients availability, microbial activity and plant growth [47].At the other end of the drainage gradient, jack pine can grow on very xeric sandy sites that are also less favourable to hardwoods, but on which ericaceous shrubs like sheep laurel can also thrive [19,48].
Both hardwood categories shared many influencing factors, which made it difficult to separate them clearly.They were more associated with well-drained sites located in warmer ecological regions located to the south of the study area.They were also less common in black spruce or jack pine ecosystems.Drainage and soil texture effects were closely linked, which explains that shade-intolerant hardwoods were more likely to be found on tills 50-100 cm thick and less on disintegration moraines and ice-contact deposits.Outwash generally feature low nutrient availability [49], which can explain why this deposit had a negative effect on the probabilities to encounter hardwoods.Similarly, previous observations have shown that paper birch, trembling aspen, pin cherry, mountain maple and speckled alder are seldom found on fluvioglacial deposits such as ice-contact deposits and outwash [50].
The link between the pre-harvest canopy abundance of intolerant hardwoods and their regeneration is not clear since some positive relationships between the abundance of paper birch in the canopy and regeneration of intolerant species were found in the RDA while negative relationships were observed in the linear mixed models.Even though paper birch can invade cutovers, it mostly regenerates through seeds that can disperse over relatively long distances [51], making it less tightly linked to the harvest site.A positive relationship was found between the presence of aspen in the canopy and regeneration of intolerant hardwoods in the linear mixed models, but not in the RDA.Since trembling aspen mainly propagates via root suckering [52], regeneration of this species is highly dependent on its presence within the canopy of the previous stand.
The abundance of non-commercial hardwoods in the regeneration layer was generally positively associated with pre-harvest hardwood stands.Although aspen and paper birch can have detrimental effects on conifer growth by intercepting light, their canopies can attenuate weather extremes and increase humidity for understory vegetation and regeneration [8,53].Their litter decomposes rapidly and increases nutrient availability [54,55].In addition, hardwoods let more light reach the understory than conifers [56], which could help semi shade-tolerant species such as mountain maple to persist under their cover.They could then expand after canopy removal.
The difficulty of differentiating between conditions conducive to a high abundance of either intolerant or non-commercial hardwoods could also be explained by the fact that the non-commercial hardwoods group included species with very different behaviours.Mountain maple and pin cherry are generally present on relatively well-drained sites [50].Speckled alder, on the other hand, is commonly associated with poorly drained sites [50].Analyzing competing vegetation by functional groups might have affected our ability to predict the presence of hardwoods like speckled alder on imperfectly or poorly drained sites.In addition, both groups could respond to common factors.
The relationship with some variables may, however, not be direct.Hence, the effect of harvesting season may be linked to the fact that wet sites, favourable to ericaceous shrubs [57] (e.g., bog Labrador tea and sheep laurel), are often selected for winter harvesting.It is also possible that the snow cover on the ground may help protect pre-established shade-tolerant regeneration during harvesting [58], making it harder for shade-intolerant hardwoods to repopulate those sites.Full tree harvesting improved the odds to encounter commercial shade-intolerant hardwoods after harvesting and had the opposite effect on ericaceous shrubs.Dragging the trees disturbs the soil, which increases soil surface temperatures and creates mineral seedbeds [59].Higher soil temperatures stimulates trembling aspen suckering [52] and the mixing of organic and mineral soil layers favour paper birch germination [60].We also observed that cut-to-length harvesting had a negative impact on the chances to encounter ericaceous shrubs after harvesting.This system has been reported to cause less damage to pre-established regeneration than full tree harvesting [59], which could help regeneration to overcome ericaceous competition.However, a given harvesting approach is usually used over a large area, so the effects of harvesting system and harvesting season may well reflect the ecosystems in which they are applied.Winter harvesting, for instance, is often preferred for wet sites where trafficability is lower and road construction costs are higher.
The effect of original disturbance could also be indirect, reflecting the vulnerability of different stand types to specific disturbances.Non-commercial hardwoods were more likely to be present after the harvesting of stands that originated from partial windthrow or partial harvesting.Hardwoods are less vulnerable to windthrow in comparison with softwoods, so partial windthrows would be more frequent in mixedwood stands and total windthrow would be more common in coniferous stands [61].Partial cutting would also be more likely in complex stands, such as mixedwood stands.Non-commercial hardwoods were probably established before these disturbances and the small gaps created could have allowed semi shade-tolerant species like mountain maple to maintain themselves in the stand until the clearcut [62,63] and then invade the cutovers [64].The negative effect of clearcutting and total windthrow on the probability of encountering commercial shade-intolerant hardwoods seems counter-intuitive at first, as the increase in light availability should have benefited these species.It is likely, however, that these disturbances were more common on cold or less fertile sites typical of black spruce or jack pine ecosystems, these being more favourable to ericaceous shrubs.

Temporal Impact of Competing Vegetation in Conditions of Marginal Regeneration
The different combinations of competing vegetation studied did not negatively affect the evolution of conifer stocking on sites with marginal regeneration, as this variable tended to improve in all combinations.In addition, significant increases of conifer free-to-grow stocking were observed on sites with either a low level of shade-intolerant hardwoods or a moderate level of ericaceous shrubs at the first inventory.In the first case, the higher availability of light provided better growth conditions for the regeneration.In the second case, this increase is likely linked to the characteristics of the sites supporting this specific combination.Most sites featuring low competition by ericaceous shrubs were imperfectly or poorly drained, while most sites featuring moderate competition by ericaceous shrubs had good or moderate drainage.Therefore, conifers happened to have better growing conditions on sites where moderate ericaceous shrubs competition occurred.Even though black spruce and balsam fir can grow on imperfectly drained sites, they have better growth rates on moderately drained sites [65].The general tendency for free-to-grow conifer stocking to improve, no matter the competition type or intensity, could be the result of the shade tolerance and nutrient requirement of the conifers.Almost every conifer encountered on the field were balsam fir or black spruce, and these species are respectively very shade-tolerant or with a broad spectrum [66] and can tolerate relatively poor nutrient levels [67].They can grow under competing vegetation covers for many decades, until they finally overtop hardwoods or ericaceous shrubs.It must be remembered here that this analysis did not include sites that were scarified and planted between the first and the second sampling.However, these treatments would normally be applied in stands that do not reach 40% conifer stocking.Since free-to-grow stocking is almost always greater than total stocking, the impact should remain minor.
Our results for sites with marginal regeneration also show that the different combinations of competition levels and composition had no effect on medium-term regeneration vigour and basal area.At the second inventory, seedlings showed good vigour based on live crown ratio and apical index.Average live crown ratio was always greater than 66% and average apical index was always greater than 1, conditions that have been associated with vigorous seedlings [25].The shade tolerance of the regenerating species could again explain the good conditions of seedlings at the second inventory.

Conclusions
The abundance and composition of competing vegetation were closely linked to broad ecosystem characteristics, namely climate, altitude, soil characteristics and dominant canopy species.Since these constitute the foundation of the ecological classification in use in Quebec [68], this classification forms a promising framework to plan efforts in competing vegetation management (e.g., mechanical site preparation, release treatment).Ericaceous shrubs were mainly found in less productive environments associated with black spruce and jack pine ecosystems, while the opposite was observed for intolerant hardwoods and non-commercial hardwoods.The low productivity could be related to either climate, altitude, soil texture or drainage.
The distinction between suitable environments for commercial shade-intolerant hardwoods and non-commercial hardwoods was less clear since they responded similarly to many variables.It is also difficult to make general statements regarding sites invaded by non-commercial hardwoods, since the diverse species included in this group represent a wide array of site requirements and regeneration strategies.
We detected some effects of harvesting system and season.However, one should be cautious in drawing general conclusions for these variables as there was some level of association between these factors and broad ecosystem characteristics in which they are preferentially applied.
Our results help to anticipate the abundance and nature of competing vegetation.One could then optimize regeneration inventories by placing less emphasis on conditions where competition levels are expected to be either very low, then requiring no tending, or very high, where tending would likely be needed.Our results also show that the impact of competition may be less critical than expected, as seedlings tended to overcome competition at marginal levels of free-to-grow stockings.
The database that served as the starting point of our study represented both a strength and a weakness.On the one hand, it provided information describing competition problems over a large section of the boreal forest of eastern Canada.On the other hand, the data were collected in an operational context, which resulted in limited details regarding non-commercial species.Non-commercial hardwoods and ericaceous shrubs were grouped into two large groups without mention of species.It would have been interesting to differentiate species since these have different ecological requirements and impact on regeneration.

Figure 1 .
Figure 1.Location of the study area and sampling sites in Quebec (Canada).Gray dots represent transects from the 2007-2011 inventory and red dots represent transects that were re-sampled in 2015.Ecological regions are those defined by [16].

Figure 1 .
Figure 1.Location of the study area and sampling sites in Quebec (Canada).Gray dots represent transects from the 2007-2011 inventory and red dots represent transects that were re-sampled in 2015.Ecological regions are those defined by[16].

Figure 2 .
Figure 2. (a) Configuration of the micro-plots in a transect used in the 2007-2011 monitoring of post-harvest regenerating stands.(b) Configuration of the plot and the microplots in a transect for the 2015 stratified inventory.

Figure 2 .
Figure 2. (a) Configuration of the micro-plots in a transect used in the 2007-2011 monitoring of post-harvest regenerating stands.(b) Configuration of the plot and the microplots in a transect for the 2015 stratified inventory.

Figure 3 .
Figure3.Redundancy analysis (RDA) ordination biplot showing the correlation between competing species groups (black arrows) and explanatory variables selected using a forward selection approach (grey arrows).Only variables with the highest contribution to axes RDA1 or RDA2 are labeled (i.e., coordinate on one axis was >90th quantile or <10th quantile of the distribution of variables coordinates on the same axis).Refer to Table2for variables description.

Figure 3 .
Figure3.Redundancy analysis (RDA) ordination biplot showing the correlation between competing species groups (black arrows) and explanatory variables selected using a forward selection approach (grey arrows).Only variables with the highest contribution to axes RDA1 or RDA2 are labeled (i.e., coordinate on one axis was >90th quantile or <10th quantile of the distribution of variables coordinates on the same axis).Refer to Table2for variables description.

Figure 4 .
Figure 4. Free-to-grow conifer stocking per strata of competing cover and over two time periods.LIH: low competing level of commercial shade-intolerant hardwoods, MIH: moderate competing level of commercial shade-intolerant hardwoods, HIH: high competing level of commercial shadeintolerant hardwoods, LES: low competing level of ericaceous shrubs, MES: moderate competing level of ericaceous shrubs, HES: high competing level of ericaceous shrubs.Bold indicate significance at α = 0.05.

Figure 4 .
Figure 4. Free-to-grow conifer stocking per strata of competing cover and over two time periods.LIH: low competing level of commercial shade-intolerant hardwoods, MIH: moderate competing level of commercial shade-intolerant hardwoods, HIH: high competing level of commercial shade-intolerant hardwoods, LES: low competing level of ericaceous shrubs, MES: moderate competing level of ericaceous shrubs, HES: high competing level of ericaceous shrubs.Bold indicate significance at α = 0.05.

Table 1 .
Distribution of transects in the stratified inventory based on dominance by competing vegetation groups and levels of competition.

Table 2 .
Description of the independent variables identified in the redundancy analysis and later used in the linear mixed model and the cumulative link mixed models.

Table 3 .
Summary of the linear mixed models predicting commercial shade-intolerant hardwood stocking as a function of harvesting and biophysical variables.

Table 4 .
Threshold coefficients and summary of the cumulative link mixed model predicting non-commercial hardwood percent cover as a function of harvesting and biophysical variables.

Table 5 .
Threshold coefficients and summary of the cumulative link mixed model predicting ericaceous shrub percent cover as a function of harvesting and biophysical variables.