Exploring the potential of two -aged white spruce plantations for the production of sawlog volume with simulations using SORTIE-ND

 martin.beland@umoncton.ca Abstract The main objective for even-aged plantation (EAP) management of producing sawlog material has driven practices towards low initial planting densities and lower post thinning densities. For semi-shade tolerant species, the resulting stand density potentially leaves enough growing space for the introduction of a second cohort of trees in the understory, making it a two-aged plantation (TAP). TAPs could have many silvicultural benefits, especially in sensitive areas where intensive treatments associated with EAPs are incompatible with local management objectives. White spruce ( Picea glauca ) is a good candidate species for modeling TAPs because it is the most widely planted tree species in Canada and has proven tolerance to understory planting. SORTIE-ND, a single-tree spatially explicit growth model was used to explore the yield of variable density and rotation length scenarios when each white spruce cohort is introduced mid rotation, compared to traditional even-aged management. All TAP scenarios tested produced more sawlog volume and more merchantable volume than equivalent densities of EAPs. The lowest density tested, 400 stems ha -1 planted every 35 years, had the highest sawlog yields (3.23 m 3 ha -1 yr -1 ). Considering smaller size products changes the optimum TAP scenario but maintains the advantage


Introduction
Plantations play a major role in providing goods and services valued by society and are widely used worldwide (Szulecka et al. 2014). Motivations for their use include increased yields of desired species, possibility to use improved seed provenances, early production of lumber and increased productivity of harvesting operations (Savill et al. 1997;Park and Wilson 2007;Thiffault et al. 2013;West 2014). However, plantations can present drawbacks that are sources of conflict between forest users. For example, traditional even-aged plantations have reduced complexity in horizontal and vertical structures compared to natural forest stands (Bowyer 2006;Brockerhoff et al. 2008;Blaha and Matjková 2012), which raises concerns regarding their capacity to maintain biodiversity (Brockerhoff et al. 2008; Barrette et al. 2014). Also, attenuation of their visual impact in the landscape is increasingly requested in forestry public participation committees, whereas vegetation management with chemical herbicides raises significant social concerns (Anonym 2014;Wyatt et al. 2011).
There is a growing interest in forest/land management alternatives to traditional plantation forestry (Edmunds and Wollenberg 2013). This has allowed trial and implementation of approaches of silviculture based on plantation/natural stands uneven-aged management (Mizunaga et al. 2010;Griess et al. 2017), mixed species, variable thinning (Carey 2003), or conversion from even-aged to uneven-aged (Kelty et al. 2003). Outcomes of these are particularly promising when ecological and economical synergies are considered. Some of these silvicultural alternatives attempt to mimic natural or old-growth forest stand structures and attributes, but others have no equivalent in nature and constitute novel socioecosystems (Lindenmayer et al. 2015).
In New Brunswick, Canada, the lumber industry is considered the cornerstone of the forestry sector, as sawlogs provide the biggest revenue per log (Roberts 2008). White spruce (Picea glauca (Moench) Voss.) is the most common species used in plantations of improved provenances, combined with site preparation and herbicide application, enabling early production of commercial size logs (Pelletier and Pitt 2008).
In an attempt to further shorten rotations, decrease the proportion of pulp wood produced, concentrate resources on fewer elite trees and decrease harvesting costs, a recent trend is to decrease planting densities from 2,500 to 1,200 stems ha -1 (pers. comm., Andrew Willett, Forest Research Director, JD Irving, ltd, Oct. 2018) and thin once, early in stand development, at low residual densities (Pelletier et Pitt 2008). At low densities, crowns of individual trees are fully developed, but stands remain below crown closure for a few years, which suggests that site occupancy is not optimal. Given that white spruce is semi-tolerant to shade (Humbert et al. 2007), and thus, commonly used for underplanting in Canada (Paquet et al. 2006), there is a potential to test intensive production of two-aged plantations (TAPs).
TAPs offer the advantages of improved control of competing vegetation and reduced need for herbicide applications (Groot 1999;Truax et al. 2000), establishment of adequate advanced regeneration prior to final harvest (Frank et al. 2018), reduced effects of potential natural disturbances (Šafránek et al. 2018), increased sustainability and productivity of stands (Pukkala et al. 2013) and higher proportion of sawlog-sized wood harvested (Struckmann 1983). Two-storied or two-aged silviculture has not been studied extensively yet otherwise than in the context of forests established through natural regeneration for non-timber objectives (e.g. Acker et al. 1998;Grover et al. 2014;Forbes et al. 2015;Heydari et al. 2015;Kabzems et al. 2016;Jaloviar et al. 2018;Truax et al. 2018). With a few exceptions (e.g. Chan et al. 2006), previous studies on the productivity of white spruce in two-aged structured stands have mainly focused on combinations with Populus sp. or mixtures of broadleaf species (Brace and Bella 1988;Pitt et al. 2015;Kabzems et al. 2016). Modeling of two-storied stands requires incorporating the effect of shading of the overstory on the understory cohort (Tatsuhara 2001).
Our objective was to use a modeling approach to explore the potential effects of density and length of rotation on merchantable volume and sawlog volume yields of even-aged and two-aged white spruce plantations in a context of intensive wood production in New Brunswick, Canada. Our general hypothesis was that there exist combinations of densities and rotations allowing the same yielding or over-yielding with TAPs compared to EAPs.

SORTIE-ND model
We used SORTIE-ND as a growth model to compare the yield of even-aged and two-aged white spruce plantations, the later being composed of two cohorts. SORTIE-ND is a partly mechanistic, partly empirical, single-tree spatially-explicit model, which strength is the consideration of the effects (mainly shading) of neighboring trees on tree development (Canham et al. 1999(Canham et al. , 2004Coates 2010;Hu and Zhu 2008;Brown et al. 2018). Model behaviors influencing tree development are documented in a "parameter file". The parameterization of white spruce in SORTIE-ND (Coates 2010) and extensive use of SORTIE-ND to estimate light availability to white spruce in understory conditions (Astrup and Larson 2006;Beaudet et al. 2011;Bose et al. 2015) justifies the use of this model to estimate white spruce growth response in two-aged plantations. We used the parameters developed by Coates (2010) for Interior spruce (P. glauca × P. engelmannii Parry) in British Columbia, Canada, that we adapted for volumes commonly reached by white spruce in New Brunswick.
Each planted tree was located in a XY coordinate system in simulated plots with wrap-around edges in the shape of a torus in order to eliminate edge effects (Gibson 2014). We used plots 300 m × 300 m to avoid a tree being considered its own competitor. The spatial distribution of planted seedlings belonging to each cohort was set on a regular square spacing, with each succeeding cohort of two-aged plantations disposed in staggered rows from the preceding one ( Figure 1). This distribution eliminates variations in growth response due to variable spacing and ensures optimal wood production (Pollack et al. 1990;Courbet et al. 2002;Johnstone and van Thienen 2011). Seedlings were assigned initial height randomly drawn from a normal distribution with mean of 43.5 cm and standard deviation of 14.5 cm (MRNF 2011), while excluding seedlings below 29.8 cm and above 57.2 cm. Seedling basal diameter was calculated based on a height/diameter ratio of 7.343. These dimensions correspond to 310 cm 3 containerized large-stock seedlings expected appropriate for understory establishment without herbicides (Thiffault et al. 2012). In order to make volume calculations correspond to planted white spruce in New Brunswick, "Maximum height" (28 m) and the curve's "Asymptotic height" (0.035) were modified from allometry parameters for Interior spruce. These new values were taken from height-diameter relationships included in the Acadian variant of the OSM/FORUS model (Hennigar et al. 2017). Figure 1. Three-dimension rendering of a two-aged white spruce plantation with a density of 400 stems ha -1 by cohort. Tree sizes (height, height to base of live crown, crown width and stem diameter at breast height) are to scale with stand averages simulated with SORTIE-ND before harvesting of a mature cohort at 70-yr rotation age.
In SORTIE-ND, the competition effect of neighbors on a given tree growth is first determined by its potential maximal diameter growth at breast-height (1.3 m). For each tree, this maximum (1.09 cm yr -1 ) is affected by model coefficients, the main being the Shading Effect Coefficient. Shading for a target tree is computed from species-specific light interception coefficients and cumulated along the sun's path during the growing season over the crowns of each neighboring trees. The characteristics of the tree stem is used to model the crown as a cylinder (Sattler and LeMay 2011). Crown size is estimated from linear relationships with tree height (Canham et al. 1999).
For two-aged stand simulations, harvesting was defined as cutting of all trees of the oldest cohort (the first cohort) at the end of the prescribed rotation. Since the model does not assign an age to trees and does not allow to discriminate trees based on the time of planting, commands given to the model set a range of diameters at breast height (dbh) corresponding to all trees belonging to the first cohort. No mortality on residual trees of the second cohort due to damage from harvesting operations was introduced in the model, as our objective was to focus on theoretical over-yielding. We set the probability of mortality for seedlings, saplings and adult trees to 0.0025, 0.0025 and 0.0005, respectively. The volume of each tree was defined in the "Tree Volume Calculator" of SORTIE-ND by integrating, using 0.5 cm increments, bole sections from the stump (height = 0.15 m) up to the top end diameter inside bark of 8.6 cm and 24.6 cm for gross merchantable volume (GMV) and sawlog volume, respectively.

Scenarios
We simulated 18 scenarios representing all combinations of 3 factors: plantation density (800, 1200 and 1600 stems ha -1 ); rotation length (50, 60 and 70 years); and plantation type (even-aged or two-aged). Densities and rotations were selected to range from the density and length of rotation in current-practice for the even-aged white spruce plantations in New Brunswick, towards lower densities and longer rotations to account for the impact of competition from the first cohorts on the growth of the secondary cohorts. For two-aged stands, planting was simulated as being performed at year 0 of the scenario (with half of total density) and introducing a second cohort (with another half of total density) within a delay of 50% of rotation length (e.g. a new cohort planted every 35 years for a 70-year rotation). Each time a new cohort was introduced in two-aged plantations, we stopped the simulation, introduced a new tree map, adjusted the lower bound of diameter range for harvest and extracted results.
We assessed treatment effects on stand yield using mean annual increment (MAI) in volume (m 3 ha -1 yr -1 ), obtained by dividing the growing stock by the age since stand origin. The MAI was calculated for GMV and sawlog volume as produced by SORTIE-ND's Tree Volume Calculator. For two-aged plantations, we used the volumes of 2 previous successive cutting cycles, once the scenarios were stabilized, divided by the length of the rotation, as a basis for calculations of MAI.

Results
When estimated for rotations between 50 and 70 years, sawlog volume and sawlog volume yield in even-aged plantations with 800 stems ha -1 ranged from 31-195 m 3 ha -1 and from 1-3 m 3 ha -1 yr -1 , respectively (Figure 2a and 2b), with production increasing with the length of the rotation. Even-aged plantations with higher densities and shorter rotations produced smaller yields. The lowest density (800 stems ha -1 ) accelerated by 6 to 15 years the recruitment of trees in the sawlog-size diameter class, compared to the higher simulated densities (1200 and 1600 stems ha -1 ).
Simulations of sawlog volume in TAPs took many rotations before stabilizing. Figure 3 exemplifies the stabilization pattern of a two-aged plantation established at a density of 800 stems ha -1 over a 70-year rotation. In this example, it took 280 years for the sawlog volumes of two successive harvests to vary by less than 5%. Hence, the following results focus on yields obtained after stabilization was reached.
For TAPs with densities varying between 800-1600 stems ha -1 and rotations between 50-70 years, yield of sawlog volume and gross merchantable volume spanned from 0-3 and 5-7 m 3 ha -1 yr -1 , respectively. Except for a scenario comprising 1600 stems ha -1 and a 50-year rotation that produced no sawlog volume, all tested scenarios of TAP produced a yield of sawlog volume as well as gross merchantable volume higher than that of even-aged plantations of equivalent density and rotation. The yield advantage of TAPs over EAPs increased with rotation length. As for EAPs, TAPs' yields of sawlog volume increased with lower densities and longer rotations (Figure 4). Gross merchantable volume responded differently, with higher yields obtained at the intermediate 1200 stems ha -1 planting density and increasing with rotation length (Figure 5).

Discussion
Our simulation results confirm the hypothesis that TAPs can provide sawlog volume yields not only comparable, but superior, to EAPs. Stand structural diversity has beneficial effects on stand growth as a result of the complementarity of resource use by the different vertical layers of vegetation having unique sets of habitat requirements (Lei et al. 2009). In our simulations, the longer the rotation, the better the sawlog volume yields for both EAP and TAP, but more so for TAPs. This is in accordance with previous work that suggested that delaying final harvest in EAPs of conifer produces higher sawlog yields (Walker 1993;MacDonald and Hubert 2002;Felton et al. 2017). Inversely, pulp volume decreases as rotation increases, be the stand even-aged or two-aged (results not shown). For similar combinations of density and rotation, we observed this trend for all comparisons between EAPs and TAPs. This was expected for TAPs; as rotation extends, the understory cohort increases its contribution to sawlog production and all stems of the new cohort become merchantable as the oldest cohort is harvested. We used 24.6 cm as a lower limit of sawlog size material at the small end of the trees for comparing scenarios on a relative basis. This size corresponds to the criteria used on NB public lands for stumpage calculations. It does not, however, factor in cull, saw kerf and other quality criteria qualifying a log for lumber. This criterion affects the time required to reach sawlog size. Not considering this had the effect of slightly overestimating yields. But more importantly, using a sawlog size threshold higher than used in other publications (e.g. Routa et al. 2019) would have delayed the production of sawlog volumes and make our estimations lower. Nevertheless, these discrepancies are affecting the sawlog volume outcomes but not the actual growth of trees in the model. Within a range of rotations far from stand breakup, which could happen as late as 240 years (Viereck et al. 1983), they should not affect the relative comparisons between scenarios.
In our simulations, sawlog volume increased as densities decreased. The densities we tested are smaller than those usually used in EAPs. Achieving the foreseen yields in sawlog volume assumes these low densities will not pose problems. Possible concerns related to low planting densities include excessive branch development and consequent lower wood quality, development of competing ground cover vegetation, and risk of windthrow. Although it remains to be studied empirically, the two-aged structure may alleviate some of these concerns since seedlings are growing in the understory of a primary cohort for the first half of the rotation. Low plantation densities (Wilson and Oliver 2000;Stark et al. 2013) and the cover of a precedent cohort can have additional dynamic effects on live crown and H/D ratios. Our simulations with SORTIE-ND however did not take these into account since the only process affected by the light environment in the model is diameter growth; height and crown architecture are derived from diameter. A more dynamic modelling of height and crown size would probably have generated larger crowns, hastening crown closure, reducing seedling growth and overyielding compared to EAPs.
Overstory cover and low planting density ought to have contradictory effects on crown development and taper. There may be a combination of planting density and rotation lengths where the use of TAPs would provide optimum light conditions that offer a compromise between growth and quality. One of the challenges of managing plantations for higher structural diversity is to maintain leaf area of the dominant layers below a threshold that still allows adequate seedlings and sapling growth in the understory (Jobidon 2000;Astrup and Larson 2006). Growth rates are mainly determined by social position within the stand (Mailly 2014). Crown development, diameter growth, and to a lesser extent, height growth [of conifers] are tightly linked to light availability (Parent and Messier 1996). As a result, lower social class trees get non-negligible impact on their branch development, especially at a young age, from shade casted by higher class neighbour crowns (Mailly 2014). Conversely, the yearly growth ring width of a dominant tree may remain constant below crown base, which tends to decrease stem taper (Newnham 1990).
In our simulations, a delay was necessary for the harvested volumes to stabilise from one cohort to the next. While the first cohort was free to grow for the whole rotation, these conditions never reoccurred as trees were introduced midrotation and were shaded by earlier cohorts in the upper strata. The fast growth rate of the first cohort allowed it to reach bigger sizes, which in turn reduced light availability in the understory to the point at which mortality occurred in the following cohorts. This phenomenon faded with successive rotations until growing conditions were stabilized between cohorts. This suggests that if TAPs were implemented, attenuation measures should be taken to hasten stabilization (e.g. reducing density of the first cohort below the tested densities of 400, 600 and 800 stems ha -1 , with thinning).
Optimal light levels tend to favor an accelerated growth in shade-tolerant conifer, while impeding the excessive recruitment of shade-intolerant species (Brace and Bella 1988;Groot 1999). By fine-tuning the light conditions in the understory and by introducing a new cohort of trees able to capture it, there is a possibility that competing vegetation be less aggressive, thus reducing the need for herbicide application. This could be a major advantage in contexts where the latter are not an option or are incompatible with management objectives. To achieve that, the understory cohort might have to be planted at a higher density to achieve crown closure more rapidly (e.g. 800 stems ha -1 instead of 400) and later be thinned concurrently with harvesting of the overstory cohort. The extra seedlings, although representing extra costs, could offer flexibility to select the best trees and to account for some of them being subject to damage due to harvesting operations.
We based growth and mortality modelling parameters on Interior spruce growing in British Columbia and adapted volumes to white spruce growing in New Brunswick. The absolute values we obtained should thus not be considered as representative of the yields that can be expected from these species in these specific regions with a two-aged plantation approach. Still, they enable relative comparisons between theoretical scenarios representing a gradient of initial stand density. However, simulated growth of EAP in the present study are similar to those reported in yield tables for white spruce plantations in Québec (Popovich 1977;Pothier and Savard 1998) and slightly lower than those reported in a thinning trial in New Brunswick (Pelletier and Pitt 2008). Other combinations of initial planting densities and rotation lengths should be tested through simulations and may allow to determine an optimum yield as in the concept of biological rotation for even-aged stands. However, this optimum might not be the best scenario for implementing TAPs in real life situations. Chosen densities and plantation patterns may, among other factors, influence wood quality (e.g. size and abundance of knots, proportion of juvenile wood). Ground vegetation, natural regeneration of trees and microsite quality may interfere with seedling growth and react differently depending on the scenario. Some situations will require site preparation or vegetation management, not accounted for by SORTIE-ND given the set of processes and parameters we used. These could entail extra costs that should be integrated into a formal financial analysis taking into account the net value of the whole assortment of wood products, actualized and repeated at perpetuity (Chang 1998). Moreover, the implementation of TAPs raises other questions: what are the actual light levels encountered in the understory of post thinning plantations? Should a denser understory cohort be planted to compensate for mortality due to damage from harvest operations or for inadequate light conditions for the survival and growth of some of the seedlings? The challenge is to reach, within the smallest delay possible, an optimum scenario maximizing sawlog yields. Ultimately, field trials of understory planting in post thinning even-aged plantations could be performed to adjust the prescription and measure other interactions not accounted for by the simulations.

Acknowledgments
The study was performed with support from Université de Moncton during the undergraduate thesis by VSG, following the original idea provided by MB. Draft revision and expertise in plantation forestry were provided by NT. Charles Canham and Lora Murphy from the Cary Institute of Ecosystems Studies provided Technical help with the use of SORTIE-ND. Detailed instructions on the use of SORTIE-ND and base parameter file for interior spruce in British Columbia (BC), Canada were provided by Dave Coates, retired research silviculturist for the BC government.