Long-term effects of harvest on boreal forest soils in relation to a remote sensing-based soil moisture index

Forest Ecology and Management 462 (2020) 117986

Contents lists available at ScienceDirect

Forest Ecology and Management journal homepage: www.elsevier.com/locate/foreco

Long-term effects of harvest on boreal forest soils in relation to a remote sensing-based soil moisture index

T

⁎

Paul D. Sewell, Sylvie A. Quideau , Miles Dyck, Ellen Macdonald Department of Renewable Resources, University of Alberta, Edmonton, AB T6G 2H1, Canada

A R T I C LE I N FO

A B S T R A C T

Keywords: Boreal mixedwood forest Soil carbon Long-term soil response Tree species effects Variable retention harvest Depth-to-water index

Storing a significant portion of the global carbon (C) stocks, soils of the boreal forest display a high degree of spatial heterogeneity across the landscape, arising from variation in forest structure and landscape morphology, as well as natural and anthropogenic disturbances. Because of this high degree of variability, accurately quantifying C storage in this ecosystem poses a challenge. Forestry is an important feature of Canada’s natural resource-based economy, but there is still considerable uncertainty on how management practices will affect boreal soil C sequestration in the long-term. With increasing pressures due to a changing climate and intensified forest management, developing better tools and techniques to quantify soil C dynamics is of paramount importance. In this study, we measured soil C stocks and associated properties in the forest floors and mineral soils (0–7 cm) of conifer-dominated and deciduous (broadleaf)-dominated boreal forest stands, 17 years following variable retention harvest. We investigated if the Wet Areas Mapping-based depth-to-water (DTW) index, derived from remotely sensed Light Detection and Ranging (LiDAR) data, was related to soil properties of unharvested stands. In addition, relationships between harvest-induced changes in soil properties and the DTW index were examined. In unharvested stands, several forest floor properties were related to DTW. However, these relationships were stand-specific. In unharvested conifer-dominated stands, forest floor C stocks were positively related to site wetness, while in harvested stands there was no relationship; this suggests that C was differentially lost from wetter sites following harvesting. Conversely, in unharvested deciduous-dominated stands, there was no relationship between forest floor C stocks and site wetness, but in harvested stands, wetter sites had higher C stocks. The DTW index was more strongly related to soil properties in the mineral soil (0–7 cm) than in the forest floor. In both forest cover types, mineral soil C and N concentrations, and C stocks, increased with increasing wetness. Relationships between mineral soil properties and DTW were not stand-specific, and were of similar magnitude under deciduous-dominated and conifer-dominated cover. In addition, harvest and forest regeneration had limited impact on the relationships between DTW and mineral soil properties. This study highlights the potential of using remote sensing and the DTW index to model forest floor and mineral soil properties in the boreal forest. In turn, this approach may be utilized in effective forest management strategies that aim to conserve boreal C stocks.

1. Introduction The boreal forest is one of the largest forest ecosystems in the world, ranging from 50 °N to the circumpolar zone (Lorenz and Lal, 2010). It contains an estimated 367–1716 Pg of carbon (C), which represents about a third of the global terrestrial C stocks (Bradshaw and Warkentin, 2015; IPCC, 2001). The boreal forest comprises a large portion of Canada, spanning 552 million ha (Brandt, 2009); indeed, this is the largest forest biome in the country. Soils are important C reservoirs in the boreal, comprising up to 85% of its total stocks (Dixon et al., 1994). In addition to its size and vast C stores, this ecosystem is

⁎

also important to Canada’s natural resource-based economy, both for the production of wood products and for oil and gas exploration (Burton et al., 2003). Combined, periodic natural and anthropogenic disturbances (e.g.; fire, insects, silviculture, climate change) have shaped the current structure of the boreal forest, both in terms of its aboveground stand dynamics, and its belowground C stores. As the extent and severity of such disturbances continues to grow, so will the importance of understanding and managing their effects on boreal C stocks. Climatic conditions in boreal forests inhibit the degradation of plant litter, allowing it to accumulate as surface organic matter. This is a

Corresponding author at: Department of Renewable Resources, 3-40 Earth Sciences Building, Edmonton, AB T6G 2E3, Canada. E-mail address: [email protected] (S.A. Quideau).

https://doi.org/10.1016/j.foreco.2020.117986 Received 5 November 2019; Received in revised form 10 February 2020; Accepted 10 February 2020 0378-1127/ © 2020 Published by Elsevier B.V.

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

storage and sequestration in the boreal, but this has yet to be demonstrated. The EMEND project is a unique experiment offering the opportunity to observe the long-term response of soil properties to disturbance and regeneration. In this study we looked at soil properties in the forest floor and upper mineral soil (0–7 cm) in conifer-dominated and deciduous-dominated stands across the EMEND project, at varying degrees of retention harvesting. Our sampling was stratified along a moisture gradient modeled using the DTW index as per Bartels et al. (2019). Our objectives were to investigate if the DTW index was related to soil properties and C stocks, and to compare these relationships between unharvested conifer-dominated and deciduous-dominated stands. We were also interested to see if variable retention harvest, and the subsequent regeneration of trembling aspen in both stand types, altered the relationships between these soil properties and the DTW index.

feature unique to boreal forest soils and underlies their significance to the global C cycle (DeLuca and Boisvenue, 2012). Across the range of the boreal there is a great deal of spatial heterogeneity arising from differences in landscape morphology and forest structure, which ultimately impacts the quantity and chemical quality of C inputs, as well as the physical constraints on decomposition or organic matter accumulation. Stands of coniferous and deciduous (i.e., broadleaf) trees have vastly different soil properties arising from differences in tree morphology and the physical and chemical properties of litter inputs into the soil (Hannam et al., 2004; Laganière et al., 2017, 2013; Lindo and Visser, 2003). In addition to biological influences, hydrological conditions play an important role on C dynamics of the boreal forest. The depth to the water table and soil moisture regime are both well-known regulators of soil C. High soil organic C content is associated with soils with poor drainage, with periodic saturation leading to anoxic conditions that inhibit decomposition processes (Rapalee et al., 1998). Additionally, dissolved organic C, which may constitute an important fraction of C stocks in boreal mineral soils (Norris et al., 2011), is thought to accumulate in low lying areas (Dalsgaard et al., 2016). In upland forest floors, organic C stocks have been found to increase not only due to limited decomposition, but also because of greater litter inputs from enhanced biomass production (Olsson et al., 2009). Interestingly, biological and hydrological constraints on C cycling do not occur independently in the boreal, as stands of increasing conifer cover generally occupy wetter sites across the boreal mixedwood forest when compared to deciduous-dominated stands (Nijland et al., 2015). Forest management has evolved over the past few decades from solely optimizing wood production to viewing forests as complex ecosystems, where management aims to maintain a broad range of ecological goods and services (Burton et al. 2003). Retention harvesting was introduced over 25 years ago, to offset the negative impacts to biodiversity brought on through clear-cutting (Lindenmayer et al., 2012). Retention harvesting involves retention of live trees during harvest in order to maintain aspects of forest structure and function and to more closely mimic natural disturbances (Fedrowitz et al., 2014). In 1999, this strategy was applied in the boreal mixedwood forest of Alberta, at the Ecosystem Management Emulating Natural Disturbance (EMEND) project. At EMEND, harvest was conducted at varying retention levels among stands of the dominant tree species across the landscape, with the ultimate objective of monitoring the uninterrupted evolution in vegetation and soil characteristics for the full lifespan of the forest (Spence et al., 1999). Logging operations at EMEND did well to preserve the soil through dedicated machine corridors and winter harvest, and previous reports found that changes in soil properties were limited six years following variable retention harvest (Kishchuk et al., 2014). In the years following harvest, regeneration was dominated by trembling aspen (Populus tremuloides Michx.) in all stand types (Echiverri, 2017; Nijland et al., 2015). Increasing C storage and sequestration in the boreal has become a management objective of increasing priority due to heightened C emissions (IPCC, 2007). Due to the complexity of C dynamics in boreal forest soils it is paramount to develop better tools and strategies to obtain more reliable information for decision making (Townshend et al., 1991). Remote sensing technologies have advanced our ability to map both topographic and vegetative features of the boreal forest (Nijland et al., 2015). The Wet Areas Mapping- (WAM) based Depth-toWater (DTW) index has been used by industry for operations planning; prior research at the EMEND project, as well as in other forested areas in Alberta, has demonstrated its utility in predicting site index, stand type, and in modeling bryophyte and understory community composition and abundance (Bartels et al., 2018; Bartels et al., 2019; Echiverri and Macdonald, 2019; Murphy et al., 2007; Oltean et al., 2016). In addition, it has been shown to outperform more conventional modelling methods when mapping soil type and forest floor depth (Murphy et al., 2011). The DTW index may also be a useful tool to model carbon

2. Materials and methods 2.1.

STUDY SITES

The study was conducted at the EMEND Project (56° 46′ 13″ N, 118° 22′ 28″ W), located in the Lower Boreal Highlands subregion of the boreal mixedwood in Alberta, Canada. EMEND is an experimental forest where green tree variable retention (VR) harvesting was conducted in 1999; for a full description of the experimental design see Spence et al. (1999). The mean annual temperature of this subregion is −1 °C, with mean annual precipitation of 495 mm, of which a third accumulates as snowfall (Natural Subregions Committee, 2006). A detailed description of the soils at EMEND can be found in Kishchuk (2004). In brief, soil development primarily occurred on glaciolacustrine and glacial till deposits. Soils are predominately fine-textured and contain very few coarse fragments. They belong to the Luvisolic and Brunisolic orders, with Gleysols and Organic soils occurring to a smaller extent in depression or discharge areas. Our study was based in conifer-dominated stands (CDOM; consisting of > 70% conifer canopy cover) primarily composed of white spruce (Picea Glauca (Moench) Voss), and deciduous-dominated stands (DDOM; consisting of > 70% broadleaf cover) of mostly trembling aspen (Populus tremuloides Michx.), and balsam poplar (Populus balsamifera L.). We sampled across various retention levels, including clear cut (2%), 20%, 50%, and lastly 100% retention, which served as an uncut control. Each retention treatment covered 10 ha in size and is herein referred to as a “compartment”. Each compartment was replicated in triplicate across the experimental area, which spans in total greater than 7800 ha. The EMEND experimental site encompasses productive sub-mesic to sub-hygric site types (Nijland et al., 2015). A moisture gradient was modelled across the experiment using the WAM-based DTW index generated from a high resolution (1 m2) digital elevation model based on LiDAR data collected in 2008 (Nijland et al., 2015). DTW serves as a measure of the probability of soil to be saturated. Its units are in meters and it approximates the depth to the water table, with lower values predicting wetter soils, and higher values predicting drier soils (Murphy et al., 2007). DTW values are sensitive to the catchment area needed to form a flow channel, known as the flow-initiation threshold, which ranged from 0.5 ha−1 to 16 ha−1 in this study. For example, DTW based on a 0.5 ha−1 flow initiation threshold will predict that a greater proportion of the landscape is wet compared to a higher threshold (Bartels et al., 2018). Sampling occurred within three replicate compartments of the four harvesting treatments (2%, 20%, 50%, and 100% retention) within each of the two stand types, conifer-dominated and deciduousdominated, for a total of 24 compartments. Within each compartment, 7 – 13 sampling locations were established along the gradient of wetness, as indicated by the DTW index as per Echiverri and Macdonald (2019). Previous work at the EMEND project reported that deciduous2

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

stands only, and model (2) was used when all retention levels were included:

dominated stands covered a wider range of DTW values, but tended to occur on the drier end of the moisture gradient, compared to coniferdominated stands (Nijland et al., 2015); in our study, DTW values at a 4 ha−1 initiation threshold ranged from 0.01 m to 11.93 m in deciduous-dominated stands (with a median value of 0.73 m), and from 0.00 m to 6.18 m in conifer-dominated stands, with a median value of 0.38 m.

Y = μ + DTW + e

(1)

and

Y = μ + DTW + Retention + DTW ′Retention + Compartment + e

(2)

where: Y is soil property, µ is the intercept, DTW is the depth-to-water index, and Retention corresponds to the harvesting retention level (2%, 20%, 50%, and 100%); these were included as fixed effects. Because of the nested structure of the sampling design (sampling locations nested within compartments), compartment was included as a random (blocking) effect; e is the random error addressing the variability between each replicate compartment. The DTW index was natural logtransformed for all analyses to reduce positive skew. A weighted variance structure, varIdent, was included to account for heterogeneity of model residuals among compartments. When the interaction between DTW and retention level was significant (α = 0.10), post-hoc analyses were conducted to compare the slope coefficients (between soil parameters and the DTW index) among retention treatments using the lsmeans package (Lenth, 2016). The conditional R2 for each model was generated using the MuMIn package (Barton, 2018). Models were developed based on the LiDAR-derived depth to water index, using DTW values computed from a 4 ha−1 flow-initiation threshold, which was selected after initial analyses in which we used the second order Akaike Information Criterion to compare among models based on seven different flow-initiation thresholds (0.5 ha−1, 1 ha−1, 2 ha−1, 4 ha−1, 8 ha−1, 12 ha−1, and 16 ha−1) for both stand types. The ΔAICc of the 4 ha−1 threshold was lower than 2.0 most frequently in both soil layers, and was also selected for modeling bryophyte assemblages previously at the EMEND project (Bartels et al., 2018; Sewell, 2018).

2.2. Soil sampling and laboratory analyses At each selected sampling location, nested plots had been previously established for research investigating the influence of DTW on understory vegetation and regrowth with variable retention harvest treatments (Bartels et al., 2018; Echiverri, 2017). In general, regrowth was higher in stands that were deciduous-dominated pre-harvest and at the drier end of the DTW gradient; this was due to rapid and abundant vegetative regeneration of aspen (Nijland et al., 2015). At each sampling location, an undisturbed and representative area was selected and cleared of live vegetation. Green moss was excluded from the samples to stay consistent with prior soil collection at EMEND (e.g.; Kishchuk et al., 2014). A plug of the entire forest floor was collected using a 100 cm2 sampling frame. A known volume of the top 0–7 cm of the mineral soil was collected using a metal core (7.3 cm internal diameter) for bulk density and chemical analysis. Samples were stored in coolers on ice until the end of the day, and then they were weighed for field moisture content and left to air dry until further laboratory analysis. Gravimetric water content was determined on both forest floor and mineral soil samples by oven drying forest floor at 65 °C for 48 h, and mineral soil at 105 °C for 24 h (Carter and Gregorich, 2008). Bulk density was calculated as the quotient of oven-dry soil weight and volume of the soil material. Oven-dried forest floor was sieved to < 4 mm and mineral soil to < 2 mm, and the fine fractions were retained for further chemical analysis. The pH of forest floor was measured from a suspension of 5 mL of soil to 25 mL of 0.01 M CaCl2 (ISO, 2005). Subsamples were ground with a ball mill at 30 Hz for 30 s (Retsch MM200). Total carbon (TC) and total nitrogen (TN) concentrations as well as the natural abundance of 13C and 15N isotopes were measured with a ThermoScientific Flash 2000 Organic Elemental Analyzer, coupled to a Delta V Advantage IRMS. The δ 13C values (‰) were referenced to the VPBD standard, and δ 15N was referenced to air.

3. Results and discussion 3.1. Forest floor relationships with DTW and the effects of harvest 3.1.1. Unharvested controls In unharvested conifer-dominated stands, forest floor thickness, C stocks, and C:N ratios increased with increased wetness (i.e.; lower DTW values), as evidenced by negative DTW slope coefficients (Table 1). In a study relating hydrology to forest soil properties in

2.3. Statistical analyses

Table 1 Outcome of linear mixed-effect models relating forest floor properties in unharvested (control) stands with the natural log-transformed Depth-to-Water (DTW) index. DTW values were calculated using the 4 ha−1 flow initiation threshold. Significant relationships (p < 0.10) are denoted by an asterisk. CDOM: conifer-dominated; DDOM: deciduous (broadleaf) dominated stands.

Statistical analyses were performed on the following forest floor and mineral soil properties: forest floor thickness, pH, total C and N concentrations, C/N ratios, C stock (Mg ha−1), and C and N stable isotope values (Tables S1 and S2) using R statistical computing software (R Core Team, 2018). To investigate the relationships between soil properties and the DTW index, we used linear mixed-effect models with the nlme package in the R statistical environment (Pinheiro et al., 2015). Analysis was initially conducted on unharvested treatments, to evaluate the relationship between soil properties and the DTW index. We then incorporated the various retention levels, to probe the effect of harvest on these relationships. These preliminary models, which included forest type, DTW, and harvest level demonstrated that there were significant interactions between forest type and DTW and/or harvest level for some soil properties. This is similar to what had been reported for relationships between vegetation attributes and the DTW index (Bartels et al., 2019; Echiverri and Macdonald, 2019). For this reason, and also because it had already been established that the conifer- and deciduousdominated forest types were located at different positions along the DTW gradient, we chose to conduct analyses separately for the two stand types (conifer-dominated and deciduous-dominated). Consequently, the following two models were applied to each stand type, where model (1) was used to assess relationships in the unharvested

Forest floor property –

Stand type –

DTW index Slope coefficient

p-value

Thickness (cm)

CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM

−1.31 0.26 0.03 −0.09 −0.25 0.04 0.03 −0.01 −1.45 −0.11 −6.26 0.23 −0.14 −0.15 0.14 −0.20

0.0024 0.2764 0.5725 0.0011 0.2881 0.8336 0.4235 0.8388 0.0180 0.6703 < 0.0001 0.7959 0.1475 0.0105 0.1441 0.0655

pH Total Carbon (%) Total Nitrogen (%) C:N ratio −1

Carbon stock (Mg ha

3

δ

13

C (‰)

δ

15

N (‰)

)

*

*

* *

* *

R2

0.45 0.36 0.20 0.91 0.53 0.18 0.07 0.56 0.49 0.33 0.56 0.50 0.04 0.64 0.19 0.82

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

Table 2 Outcome of linear mixed-effect models relating forest floor properties to the natural log-transformed Depth-to-Water (DTW) index, retention level (clearcut, 20% or 50% retention, and unharvested control), and their interactions. DTW values were calculated using the 4 ha−1 flow initiation threshold. Significant differences in slope coefficients (response variable as a function of DTW) among retention levels are indicated by different lowercase letters (p < 0.10). CDOM: conifer dominated; DDOM: deciduous dominated. Forest floor property –

Stand Type

DTW index

Retention level

Interaction (DTW index * retention level)

–

Slope coefficient

p-value

p-value

p-value

Clearcut Slope coefficient

20% Retention Slope coefficient

50% Retention Slope coefficient

Control Slope coefficient

R2

Thickness (cm)

CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM

−0.36 −0.35 0.06 −0.09 0.36 −0.26 −0.06 −0.05 0.91 0.16 −0.39

0.0045 0.0004 0.1271 < 0.0001 0.0009 0.3983 0.2235 0.0161 0.0537 0.0308 0.0060

0.8783 0.6995 0.0981 0.9817 0.2324 0.4797 0.7596 0.6433 0.4555 0.5748 0.9786

0.1704 0.1032 0.7188 0.6599 0.0288 0.1865 0.1311 0.5884 0.0042 0.5958 0.0138

−0.36 −0.35 0.06 −0.09 0.36 ab −0.26 −0.06 −0.05 0.91 a 0.16 −0.39 a

−0.46 −0.42 0.03 −0.06 0.55 ab −0.05 0.00 −0.01 0.01 ab 0.09 −1.25 ab

−0.16 −0.28 0.00 −0.04 0.64 a 0.33 −0.03 −0.03 0.62 a 0.24 −0.55 a

−1.19 0.24 0.03 −0.09 −0.23b 0.06 0.04 −0.01 −1.45b −0.11 −6.14b

0.06 0.74 0.05 0.73 0.51 0.09 0.35 0.68 0.53 0.63 0.31

DDOM CDOM DDOM CDOM DDOM

−2.36 0.04 0.01 −0.15 0.08

< 0.0001 0.9203 0.3911 0.9341 0.3184

0.8168 0.1021 0.8434 0.2713 0.8492

0.0243 0.4206 0.0520 0.0056 0.0161

−2.36 ab 0.04 0.01 a −0.15b 0.08 ab

−3.64b −0.04 0.03 ab 0.10 a −0.02 ab

−1.52 ab 0.00 −0.01 ab −0.06 ab 0.12 a

0.01 a −0.14 −0.15b 0.15 a −0.21b

0.70 0.05 0.82 0.16 0.69

pH Total Carbon (%) Total Nitrogen (%) C:N ratio Carbon stock (Mg ha−1) δ

13

C (‰)

δ

15

N (‰)

organic matter is enriched in 13C and 15N compared to fresh litter inputs, and that isotopic enrichment increases with greater decomposition (Natelhoffer and Fry, 1988; Preston et al., 2009; Quideau et al., 2003). However, in our study, δ13C and δ15N values were most enriched at the wetter end of the DTW gradient; i.e., where decomposition was likely the slowest. Similarly, previous research demonstrated that a common response with increasing wetness is lower plant δ13C values (Farquhar et al., 1989), which again would explain a relationship with DTW opposite to what we observed. Consequently, a more likely reason for the change in isotopic abundance along the DTW gradient would be a shift in vegetation communities with distinct isotopic values. In a study on the isotopic composition of various boreal plant species, Brooks et al. (1997) reported wide variation in natural abundance 13C among species; more specifically, Flanagan et al. (1997) reported greater discrimination values (i.e.; lower δ13C values) for shrubs and forbs when compared to deciduous trees. At EMEND, the drier end of the DTW index was associated with greater abundance (cover) of understory vascular plants; in conifer-dominated stands this was largely due to increased forb cover while on deciduous-dominated stands it was attributed to increased shrub cover (Echiverri and Macdonald, 2019). From this we infer that an increase in litter inputs from these understory plants was likely driving the changes in δ13C and δ15N values along the modelled moisture gradient (Table 1).

Sweden, Olsson et al. (2009) similarly reported an increase in both forest floor thickness and C stocks with increases in soil water levels, suggesting that this accumulation was due to an increase in productivity, litter quality changes, and/or limited decomposition due to anoxic conditions during periodic soil saturation. Although overall site productivity at EMEND was found to increase towards the dry end of the DTW gradient (Nijland et al., 2015), bryophyte cover increased at wetter sites in conifer-dominated stands (Bartels et al., 2018). Bryophytes are known to constitute a significant proportion of forest floor C in conifer-dominated forests (O’Connel et al., 2003), and to have a relatively high C:N ratio compared to herbaceous plants (DeLuca and Boisvenue, 2012). In addition, C:N ratios typically decrease with decomposition in both litters and bulk forest floors (Preston et al., 2009). At EMEND, high soil water levels may have limited decomposition, leading to the higher C:N ratios observed towards the wet end of the gradient (Table 1). Further, the decomposition rate of bryophyte litter is notably low compared to other understory vegetation due to the presence of recalcitrant compounds, and compounds inhibiting microbial activity (Cornwell et al., 2008; DeLuca and Boisvenue, 2012; Harden et al., 1997). Lastly, high bryophyte cover could intensify the effects of high water on decomposition, by lowering soil temperature and oxygen levels (Oechel and Van Cleve, 1986). Consequently, in our study, the observed increase in forest floor thickness, C stocks and C/N ratios with wetness was likely related to both anoxic conditions restricting C turnover, and an increase in bryophyte abundance. In contrast to the conifer-dominated stands, carbon stocks and forest floor thickness were not related to the DTW index in the deciduousdominated forests (Table 1). Nijland et al. (2015) previously reported that deciduous-dominated stands occupied the driest sites of all forest types at EMEND. In our study, the mean DTW in deciduous-dominated stands was 2.01 m, while conifer-dominated stands occupied wetter sites, with a mean DTW of 0.83 m, although deciduous-dominated stands occupied a wider range of moisture levels. This may have given rise to competing forest processes in the deciduous-dominated stands; for example, the accumulation of forest floor due to reduced decomposition with increasing wetness may have been offset at the dry end of the gradient by the increased productivity (Nijland et al., 2015), and greater litter inputs. In unharvested deciduous-dominated stands, increasing wetness was associated with increases in forest floor pH, and δ13C and δ15N values (Table 1). Previous work in forest ecosystems has shown that soil

3.1.2. Harvested conifer-dominated stands Results from the mixed model analysis revealed significant interactions of the retention level and DTW (Table 2), indicating that harvesting, and the subsequent aspen regeneration, altered the relationships between DTW and each soil property. However, there was little difference among the retention treatments themselves. In several instances, the slope coefficient (response variable as a function of DTW) for a harvested treatment differed from the control, but there was almost no difference in the slope coefficients among the three harvesting levels (listed under Interaction in Table 2). For example, the forest floor C:N ratio showed a negative slope coefficient (-1.45) in the control treatment (Table 2), indicating that it increased towards the wetter end of the DTW gradient (Fig. 1). However, in harvested stands, the opposite was found; the C:N ratios decreased with increasing site wetness (Fig. 1, Table 2). Similarly to the C:N ratio and total C concentrations, the relationship between forest floor C stocks and DTW was altered by harvest and 4

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

aspen regeneration as evidenced by the significant interaction terms in Table 2. In the uncut controls, C stocks sharply rose at wetter sites (slope coefficient of −6.14, Table 2), while in the harvested stands there was only a very weak increase with wetness, with slope coefficients ranging between −0.39 and −1.25 (Fig. 2, Table 2). This change reflects losses of forest floor C at the wet end of the gradient in harvested stands, likely due to enhanced decomposition through disruption of the canopy and increases in soil temperatures (Lloyd and Taylor, 1994). While bryophyte cover increased towards the wet end of the gradient under all treatments, harvest led to significant losses in cover overall (Bartels et al., 2019). In addition, decreases in both C concentrations and C:N ratios towards wetter sites in harvested stands (Table 2, Fig. 1), further indicates enhanced decomposition (Preston et al., 2009). On the other end, at drier sites in harvested stands, there was an increase in C stocks relative to the undisturbed controls (Fig. 2), which likely resulted from the enhanced aspen regeneration towards the dry end of the gradient at EMEND (Nijland et al., 2015). 3.1.3. Harvested deciduous-dominated stands Similar to conifer-dominated stands, relationships between DTW and forest floor properties in deciduous stands varied significantly between the uncut control and harvested treatments, but not among the retention treatments themselves (Table 2). In uncut deciduous-dominated stands, δ13C and δ15N values were negatively related to DTW, increasing at wetter sites (Fig. 3; Table 2). Harvesting and regeneration altered these relationships, which became very weak, with slopes approaching zero (Fig. 3; Table 2). Previous studies at EMEND reported that total understory cover and community composition were altered by harvest (Echiverri and Macdonald, 2019; Echiverri, 2017), which may have shifted δ13C and δ15N values. In addition, increases in δ13C and δ15N values at the dry end of the gradient following harvest may have resulted from a rise in soil temperature and enhanced decomposition (Natelhoffer and Fry, 1988). Lastly, the proliferation of trembling aspen

Fig. 1. Relationship between forest floor C:N ratios in harvested and unharvested (control) conifer-dominated stands and the natural logarithm of the Depth-to-Water (DTW) index (smaller DTW values indicate wetter sites). Regressions lines from the linear-mixed model analysis are shown for each retention harvest treatment, with dotted lines depicting harvested treatments, and the solid line corresponding to unharvested control stands. See Table 2 for a summary of significant differences in slope coefficients among harvest treatments.

Fig. 2. Relationship between forest floor carbon stocks (Mg ha−1) in harvested and unharvested (control) conifer-dominated stands and the natural-logarithm of the Depth-to-Water (DTW) index (smaller DTW values indicate wetter sites). Regressions lines from the linear-mixed model analysis are shown for each retention harvest treatment, with dotted lines depicting harvested treatments, and the solid line corresponding to control stands. See Table 2 for a summary of significant differences among slope coefficients for each harvest treatment.

Fig. 3. Relationship between forest floor δ13C values (‰) in harvested and unharvested (control) deciduous-dominated stands and the natural logarithm of the Depth-to-Water (DTW) index (smaller DTW values indicate wetter sites). Regressions lines from the linear-mixed model analysis are shown for each retention harvest treatment, with dotted lines depicting harvested treatments, and the solid line corresponding to control stands. See Table 2 for a summary of significant differences among slope coefficients for each harvest treatment. 5

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

Table 3 Outcome of linear mixed-effect models relating mineral soil properties in unharvested (control) stands with the natural log-transformed Depth-to-Water (DTW) index. DTW values were calculated using the 4 ha−1 flow initiation threshold. Significant relationships (p < 0.10) are denoted an asterisk. CDOM: conifer-dominated; DDOM: deciduous dominated. Mineral soil property –

Stand Type –

DTW index Slope coefficient

p-value

R2

Total Carbon (%)

CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM CDOM DDOM

−0.58 −0.90 −0.03 −0.09 −0.39 −0.15 −3.27 −3.08 0.006 −0.07 0.17 * −0.17

0.0106 < 0.0001 0.0923 < 0.0001 0.1948 0.5042 0.0018 0.0077 0.8755 0.1532 0.0066 0.0537

0.62 0.45 0.47 0.58 0.38 0.01 0.68 0.40 0.63 0.16 0.90 0.56

Total Nitrogen (%) C:N ratio Carbon stock (Mg ha−1) δ

13

C (‰)

δ

15

N (‰)

* * * *

* *

*

is shown by the fact that conditional R2 values were typically much higher for the mineral soil than for the forest floor models (Tables 1 and 3). In the forest floor models, there were notable differences in the slope coefficients between the two forest types (Table 1), highlighting the influence of above ground biomass on forest floor properties. In contrast, the slope coefficients for the mineral soil models were quite similar for both forest types (Table 3). Specifically, mineral soil TC and TN concentrations, and C stocks increased towards wetter sites for both forest types; i.e., as evidenced by the negative slope coefficients relating these properties to DTW (Table 3). This increase in soil C stocks with wetness is a well reported phenomenon in the boreal forest, where soil saturation generates anoxic conditions that inhibit decomposition, allowing for organic material to accumulate (DeLuca and Boisvenue, 2012; Ping et al., 2010; Rapalee et al., 1998). The fact that the two forest types were similar in terms of the relationship between mineral soil properties and DTW suggests that the influence of aboveground biomass on mineral soil properties was smaller than their link to DTW. This also concurred in Norwegian forest soils (Dalsgaard et al., 2016), where the effects of poor drainage overruled the effects of site productivity on total soil C stocks, in part due to the transport of DOC to low lying areas.

Fig. 4. Relationship between forest floor carbon stocks (Mg ha−1) in harvested and unharvested (control) deciduous-dominated stands and the natural logarithm of the Depth-to-Water (DTW) index (smaller DTW values indicate wetter sites). Regressions lines from the linear-mixed model analysis are shown for each retention harvest treatment, with dotted lines depicting harvested treatments, and the solid line corresponding to control stands. See Table 2 for a summary of significant differences among slope coefficients for each harvest treatment.

regeneration at the dry end of the gradient following harvest at EMEND (Echiverri, 2017; Nijland et al., 2015), could have reduced soil moisture levels sufficiently to induce drought stress on trees and elevate the δ13C composition of litter fall, as reported by Garten and Taylor (1992) in a temperate deciduous forest. Harvest and post-harvest regeneration influenced the relationship between DTW and forest floor C stocks in the deciduous-dominated stands. In unharvested stands there was no relationship between DTW and forest floor C stocks (with a slope coefficient of 0.01), but harvesting led to a substantial increase in stocks towards the wet end of the gradient, and stocks nearly doubled compared to the uncut control (Table 2; Fig. 4). While post-harvest regrowth of aspen at EMEND was greatest at drier sites (Nijland et al., 2015), Echiverri (2017) reported that harvest increased total understory cover at the wet end of the gradient. Hence, the increase in forest floor C stocks at the wetter sites could reflect this increase in understory cover. Alternatively, limited aspen regeneration at the wet end of the gradient (Nijland et al., 2015), may have led to a rise in the water table through reduced transpiration (Marcotte et al., 2008), leading to an accumulation of C stocks through restricted decomposition. In addition, fine root biomass typically increases with site wetness in the boreal forest, and may have led to an increase in belowground inputs, despite the limited above-ground regeneration (Kleja et al., 2008; Olsson et al., 2009; Snedden, 2013). In all cases, either alone or in combination, these processes led to the observed increase in forest floor C stocks, as well as to the decrease in δ13C and δ15N values with increased wetness in the harvested deciduousdominated stands at EMEND (Figs. 3 and 4; Table 2).

3.2.2. Harvested stands Harvest and aspen regeneration had only weak effects on mineral soil properties (few significant effects of retention or significant retention by DTW interactions, Table 4). Although relationships between DTW, and δ15N values and C:N ratios were affected by harvesting (significant DTW and retention interaction) for both forest types, the remaining soil properties were resilient to harvest. Total C, TN, and C stocks increased towards the wetter end of the DTW gradient in harvested stands, as they did in un-harvested stands, for both forest cover types (Table 4). Previous research has shown that disturbances frequent in the boreal, such as fire, harvest, and salvage logging all have longterm effects on soil properties. However, these effects are most prominent in the forest floor, and decrease with depth (Kishchuk et al., 2015, 2014; Seedre et al., 2011). 3.3. Implications for carbon management in the boreal forest Many studies on the impacts of harvesting on soil properties report significant decreases in forest floor C stocks (Clark et al., 2015). Using the DTW index, we found that losses in forest floor C in conifer-dominated stands were greatest at the wet end of the gradient (Fig. 2). Our findings suggest that these losses may be aligned with the soil moisture regime, and that in the future, the DTW index may be a useful tool to

3.2. Mineral soil relationships with DTW and the effects of harvest 3.2.1. Unharvested controls Relationships between mineral soil properties and the DTW index were stronger than those between forest floor properties and DTW. This 6

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

Table 4 Outcome of linear mixed-effect models relating mineral soil (0 – 7 cm) properties and the natural log transformed Depth-to-Water (DTW) index, retention level (clearcut, 20% or 50% retention, and unharvested control), and their interaction. DTW values were calculated using the 4 ha−1 flow initiation threshold. Significant differences in the slope coefficients (soil property as a function of the DTW index) among retention levels are indicated by different lowercase letters (p < 0.10). CDOM: conifer dominated; DDOM: deciduous dominated. Mineral soil property –

Stand type

DTW index

Retention level

Interaction (DTW index * retention level)

–

Slope coefficient

p-value

p-value

p-value

Clearcut Slope coefficient

20% Retention Slope coefficient

50% Retention Slope coefficient

Control Slope coefficient

R2

Total Carbon (%)

CDOM DDOM CDOM DDOM CDOM DDOM CDOM

−0.46 −0.30 −0.04 −0.03 0.65 −0.03 −2.22

< 0.0001 < 0.0001 < 0.0001 < 0.0001 0.5720 0.0206 < 0.0001

0.635 0.668 0.884 0.582 0.8894 0.0638 0.479

0.475 0.377 0.743 0.112 0.0150 0.0486 0.282

−0.46 −0.30 −0.04 −0.03 0.65 a −0.03 ab −2.22

−0.77 −0.45 −0.05 −0.03 −0.15 ab −0.07b −3.77

−0.55 −0.54 −0.04 −0.06 −0.11b 0.42 a −2.39

−0.60 −0.90 −0.03 −0.09 −0.40b −0.02 ab −3.15

0.56 0.31 0.47 0.64 0.60 0.23 0.36

DDOM CDOM DDOM CDOM DDOM

−1.32 −0.04 −0.05 −0.05 0.25

< 0.0001 0.0262 0.3824 0.4488 0.9567

0.445 0.404 0.754 0.959 0.802

0.664 0.458 0.395 0.005 0.002

−1.32 −0.04 −0.05 −0.05b 0.25 a

−1.75 −0.07 0.00 0.15 ab 0.11 a

−1.65 −0.03 0.02 0.02 ab −0.13b

−3.07 0.01 −0.08 0.16 a −0.17b

0.25 0.68 0.32 0.84 0.45

Total Nitrogen (%) C:N ratio Carbon stock (Mg ha−1) δ

13

C (‰)

δ

15

N (‰)

Acknowledgements

guide retention in order to better manage for C storage. Furthermore, our results indicated that in the case of the aspen-dominated stands, harvesting led to an increase in forest floor C stocks; similar to losses in the conifer-dominated stands, increases under aspen were related to the soil moisture regime, again underlining the potential for using the DTW index for forest C management. Soil relationships with the topographic DTW index were more strongly expressed in the mineral soil (0–7 cm), compared to the forest floor. Mineral soil TN, TC, and C stocks of both forest cover types increased with wetness, supporting previous research in the boreal. In contrast to results for the forest floor, harvest and aspen regeneration had little effect on relationships present in the mineral soil, highlighting its resilience to changes in aboveground biomass and forest structure brought about by harvest (Table 4). In a meta-analysis comparing stem only harvest to whole tree harvest, Wan et al. (2018) reported that while residue retention led to increases in mineral soil carbon stocks in coarse and medium textured soils, there were no differences in fine textured soils; they proposed that the increased organic matter left through residue retention led to increased turnover though a priming mechanism. The soils at EMEND are predominately fine textured, and our study similarly highlights the ability of such soils to buffer environmental changes. In terms of effective forest management, there are different needs to be met. While the focus of this work was on carbon storage, managing for habitat and biodiversity is critical for sustainable forestry. Forest floor C stocks at the wet end of the gradient in conifer-dominated stands were most sensitive to the long-term effects of harvest. Similarly, Bartels et al. (2019) reported that leaving retention at wetter sites would be most effective at maintaining bryophyte cover in this stand type, and Echiverri (2017) also suggested that retention be placed at the wet end of the gradient, as the understory vegetation in these areas was the most sensitive to the long-term effects of harvest. With these findings in mind, retention should be placed at the wet end of the gradient in conifer-dominated stands. On the other hand, in deciduous-dominated stands, the long-term effects of harvest led to an increase in forest floor C stocks at the wet end of the gradient. Echiverri (2017) recommended that retention be placed at the dry end of the gradient in deciduous-dominated stands, as these areas were less resilient to harvesting activities. In this case, the long-term effects of harvest at wetter sites may serve to increase C stores and preserve understory biodiversity.

Main support for this work came from the Natural Sciences and Engineering Research Council (NSERC) of Canada through a Strategic Grant (no. 150691074) to Ellen Macdonald, Sylvie Quideau, and Miles Dyck. It would not have been possible without support to EMEND from multiple funding sources over the years: the Sustainable Forest Management Network, Alberta Sustainable Resource Development (now Alberta Agriculture and Forestry), Daishowa-Marubeni International Ltd., Canadian Forest Products Ltd., Natural Resources Canada – Canadian Forest Service, Manning-Diversified Forest Products, Weyerhaeuser, and the Forest Resource Improvement Association of Alberta. Additional support came from a University of Alberta Northern Research Award to Paul Sewell. Lastly, Paul Sewell’s stipend was in part supported by a Queen Elizabeth II Scholarship from the University of Alberta. We gratefully acknowledge the contributions of Jae Ogilvie, Barry White, and Alberta Agriculture and Forestry for providing the lidar data, for their efforts to develop and extensively test the wet-areas mapping process in the EMEND area, and for providing the outputs from the wet-areas mapping. For assistance with field sampling we thank Maksat Igdyrov, and for their advice on mixedmodeling analysis, we thank Laureen Echiverri, and Samuel Bartels. Appendix A. Supplementary material Supplementary data to this article can be found online at https:// doi.org/10.1016/j.foreco.2020.117986. References Bartels, S.F., Caners, R.T., Ogilvie, J., White, B., Macdonald, S.E., 2018. Relating bryophyte assemblages to a remotely sensed depth-to-water index in Boreal forests. Front. Plant Sci. 9, 1–11. Bartels, S.F., James, R.S., Caners, R.T., Macdonald, S.E., 2019. Depth-to-water mediates bryophyte response to harvesting in boreal forests. J. Appl. Ecol. 56, 1256–1266. Barton, K., 2018. MuMIn: Multi-Model Inference. R Package Version 1.15.6. Available at: http://cran.r-project.org/package=MuMIn (Accessed May 19, 2018). Bradshaw, C.J.A., Warkentin, I.G., 2015. Global estimates of boreal forest carbon stocks and flux. Glob. Planet. Change 128, 24–30. Brandt, J.P., 2009. The extent of the North American boreal zone. Environ. Rev. 17, 101–161. Brooks, J.R., Flanagan, L.B., Buchmann, N., Ehleringer, J.R., 1997. Carbon isotope composition of boreal plants: functional groupings of life forms. Oecologia 110, 301–311. Burton, P.J., Messier, C., Weetman, G.F., Prepas, E.E., Adamowicz, W.L., Tittler, R., 2003. The current state of boreal forestry and the drive for change in Towards Sustainable Management of the Boreal Forest. Eds Burton, P.J., Messier, C., Smith, D.W., Adamowicz, W.L., NRC Research Press, Ottawa, Ont 1-40.

7

Forest Ecology and Management 462 (2020) 117986

P.D. Sewell, et al.

Lindenmayer, D.B., Franklin, J.F., Lõhmus, A., Baker, S.C., Bauhus, J., Beese, W., Brodie, A., Kiehl, B., Kouki, J., Pastur, G.M., Messier, C., Neyland, M., Palik, B., SverdrupThygeson, A., Volney, J., Wayne, A., Gustafsson, L., 2012. A major shift to the retention approach for forestry can help resolve some global forest sustainability issues. Conserv. Lett. 5, 421–431. Lindo, Z., Visser, S., 2003. Microbial biomass, nitrogen and phosphorus mineralization, and mesofauna in boreal conifer and deciduous forest floors following partial and clear-cut harvesting. Can. J. For. Res. 33, 1610–1620. Lloyd, J., Taylor, J.A., 1994. On the temperature dependence of soil respiration. Funct. Ecol. 8, 315–323. Lorenz, K., Lal, R., 2010. Carbon sequestration in forest ecosystems. Springer Sci. Business Media 289. Marcotte, P., Roy, V., Plamondon, A.P., Auger, I., 2008. Ten-year water table recovery after clearcutting and draining boreal forested wetlands of eastern Canada. Hydrol. Process. 22, 4163–4172. Murphy, P.N.C., Ogilvie, J., Connor, K., Arp, P.A., 2007. Mapping wetlands: a comparison of two different approaches for New Brunswick, Canada. Wetlands 27, 846–854. Murphy, P.N.C., Ogilvie, J., Meng, F., White, B., Bhatti, J.S., Arp, P.A., 2011. Modelling and mapping topographic variations in forest soils at high resolution: a case study. Ecol. Model. 222, 2314–2332. Natelhoffer, K.J., Fry, B., 1988. Controls on natural nitrogen-15 and carbon-13 abundances in forest soil organic matter. Soil Sci. Soc. Am. J. 52, 1633–1640. Natural Regions Committee, 2006. Natural Regions and Subregions of Alberta. Publication T/852. Edmonton: Government Alberta 254. Nijland, W., Coops, N.C., Macdonald, S.E., Nielsen, S.E., Bater, C.W., White, B., Ogilvie, J., Stadt, J., 2015. Remote sensing proxies of productivity and moisture predict forest stand type and recovery rate following experimental harvest. For. Ecol. Manage. 357, 239–247. Norris, C.E., Quideau, S.A., Bhatti, J.S., Wasylischen, R.E., 2011. Soil carbon stabilization in jack pine stands along the boreal forest transect case study. Glob. Change Biol. 17, 480–494. O’Connel, K.E., Gower, S.T., Norman, J.M., 2003. Net ecosystem production of two contrasting boreal black spruce communities. Ecosystems 6, 248–260. Oechel, W.C., Van Cleve, K., 1986. The role of bryophytes in nutrient cycling in the Taiga. Forest Ecosyst. Alaskan Taiga 121–137. Olsson, M., Erlandsson, M., Lundin, L., Nilsson, T., Nilsson, Å., Stendahl, J., 2009. Organic carbon stocks in Swedish Podzol soils in relation to soil hydrology and other site characteristics. Silva Fenn. 43, 209–222. Oltean, G.S., Comeau, P.G., White, B., 2016. Carbon isotope discrimination by Picea glauca and Populus tremuloides is related to the topographic depth to water index and rainfall. Can. J. For. Res. 46, 1225–1233. Ping, C.L., Michaelson, G.J., Kane, E.S., Packee, E.C., Stiles, C.A., Swanson, D.K., Zaman, N.D., 2010. Carbon stores and biogeochemical properties of soils under black spruce forest. Alaska. Soil Sci. Soc. Am. J. 74, 969–978. Pinheiro, J.C., Bates, D.M., DebRoy, S., Sarkar, D., and R Core Team. 2015. nlme: Linear and nonlinear mixed effects models. Available at: http:// CRAN.R-project.org/ package=nlme. (Accessed on April 1, 2017). Preston, C.M., Nault, J.R., Trofymow, J.A., Smyth, C., 2009. Chemical changes during 6 years of decomposition of 11 litters in some Canadian forest sites. Part 1. elemental composition, tannins, phenolics, and proximate fractions. Ecosystems 12, 1053–1077. Quideau, S.A., Graham, R.C., Feng, X., Chadwick, O.A., 2003. Natural isotopic distribution in soil surface horizons differentiated by vegetation. Soil Sci. Soc. Am. J. 67, 1544–1550. R Core Team, 2018. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. Rapalee, G., Trumbore, S.E., Davidson, E.A., Harden, J.W., Veldhuis, H., 1998. Soil carbon stocks and their rate of accumulation and loss in a boreal forest landscape. Global Biogeochem. Cycles 12, 687–701. Seedre, M., Shrestha, B.M., Chen, H.Y.H., Colombo, S., Jõgiste, K., 2011. Carbon dynamics of North American boreal forest after stand replacing wildfire and clearcut logging. J. For. Res. 16, 168–183. Sewell, P.D., 2018. Seeing the forest for the soil: topographic controls on soi carbon dynamics in the boreal mixedwood forest. M.Sc. Thesis. In: University of Alberta, Edmonton, ABpp. 105. Snedden, J., 2013. The root distribution , architecture , transpiration and root sapflow dynamics of mature trembling aspen (Populus tremuloides) growing along a hillslope. M.Sc. Thesis. University of Alberta, Edmonton, AB. 103 pp. Spence, J., Volney, W., Lieffers, V., Weber, M., Luchkow, S., Vinge, T., 1999. The sustainable forest managment network conference science and practice: sustaining the Boreal Forest Edmonton, Alberta. Sustain. For. Manag. Netw. Conf. Sci. Pract. Sustain. Boreal For. 583–591. Townshend, J., Justice, C., Li, W., Gurney, C., McManus, J., 1991. Global land cover classification by remote sensing: present capabilities and future possibilities. Remote Sens. Environ. 35, 243–255. Wan, X., Xiao, L., Vadeboncoeur, M.A., Johnson, C.E., Huang, Z., 2018. Response of mineral soil carbon storage to harvest residue retention depends on soil texture: a meta-analysis. For. Ecol. Manage. 408, 9–15.

Carter, M.R., Gregorich, E.G., 2008. Soil Sampling and Methods of Analysis, second ed. CRC Press, Boca Raton, Fl, USA, pp. 1224. Clark, N., Gundersen, P., Jönsson-Belyazid, U., Kjønaas, O.J., Persson, T., Sigurdsson, B.D., Stupak, I., Vesterdal, L., 2015. Influence of different tree-harvesting intensities on forest soil carbon stocks in boreal and northern temperate forest ecosystems. For. Ecol. Manage. 351, 9–15. Cornwell, W.K., Cornelissen, J.H.C., Amatangelo, K., Dorrepaal, E., Eviner, V.T., Godoy, O., Hobbie, S.E., Hoorens, B., Kurokawa, H., Pérez-Harguindeguy, N., Quested, H.M., Santiago, L.S., Wardle, D.A., Wright, I.J., Aerts, R., Allison, S.D., Van Bodegom, P., Brovkin, V., Chatain, A., Callaghan, T.V., Díaz, S., Garnier, E., Gurvich, D.E., Kazakou, E., Klein, J.A., Read, J., Reich, P.B., Soudzilovskaia, N.A., Vaieretti, M.V., Westoby, M., 2008. Plant species traits are the predominant control on litter decomposition rates within biomes worldwide. Ecol. Lett. 11, 1065–1071. Dalsgaard, L., Lange, H., Strand, L.T., Callesen, I., Borgen, S.K., Liski, J., Astrup, R., 2016. Underestimation of boreal forest soil carbon stocks related to soil classification and drainage. Can. J. For. Res. 46, 1413–1425. DeLuca, T.H., Boisvenue, C., 2012. Boreal forest soil carbon: distribution, function and modelling. Forestry 85, 161–184. Dixon, R.K., Brown, S., Houghton, R.A., Solomon, A.M., Trexler, M.C., Wisniewski, J., Brown, S., Houghton, R.A., Trexier, M.C., Wisniewski, J., 1994. Carbon pools and flux of global forest ecosystems. Science 263, 185–190. Echiverri, L., 2017. Identifying understory diversity and resilience patterns with the depth- to-water index in boreal mixedwood forests. M.Sc. Thesis. In: University of Alberta, Edmonton, ABpp. 91. Echiverri, L., Macdonald, S.E., 2019. Utilizing a topographic moisture index to characterize understory vegetation patterns in the boreal forest. For. Ecol. Manage. 447, 35–52. Farquhar, G.D., Ehleringer, J.R., Hubick, K.T., 1989. Carbon isotope discrimination and photosynthesis. Annu. Rev. Plant Biol. 40, 503–537. Fedrowitz, K., Koricheva, J., Baker, S.C., Lindenmayer, D.B., Palik, B., Rosenvald, R., Beese, W., Franklin, J.F., Kouki, J., Macdonald, E., Messier, C., Sverdrup-Thygeson, A., Gustafsson, L., 2014. Can retention forestry help conserve biodiversity? A metaanalysis. J. Appl. Ecol. 51, 1669–1679. Flanagan, L.B., Brooks, J.R., Ehleringer, J.R., 1997. Photosynthesis and carbon isotope discrimination in boreal forest ecosystems: a comparison of functional characteristics in plants from three mature forest types. J. Geophys. Res. 102 (D24), 28861–28869. Garten, C.T., Taylor, G.E., 1992. Foliar δ13C within a temperate deciduous forest: spatial, temporal, and species sources of variation. Oecologia 90, 1–7. Hannam, K.D., Quideau, S.A., Oh, S.-W., Kishchuk, B.E., Wasylishen, R.E., 2004. Forest floor composition in aspen- and spruce-dominated stands of the Boreal Mixedwood forest. Soil Sci. Soc. Am. J. 68, 1735–1743. Harden, J.W., O’Neill, K.P., Trumbore, S.E., Veldhuis, H., Stocks, B.J., 1997. Moss and soil contributions to the annual net carbon flux of a maturing boreal forest. J. Geophys. Res. [Atmos.] 102, 805–816. Parry, M.L., Canziani, O.F., Palutikof, J.P., van der Linden, P.J., Hanson, C.E., 2007. Climate change 2007: impacts, adaptation and vulnerability. Contribution of Working Group II to the fourth assessment report of the Intergovernmental Panel on Climate Change. Eds. Cambridge University Press, Cambridge, UK, pp. 976. IPCC, 2001. Climate Change 2001: The scientific basis. Contribution of Working Group I to the third assessment report of the Intergovernmental Panel on Climate Change. Eds Houghton, J.T., Y. Ding, D.J. Griggs, M. Noguer, van der Linden, P.J., Dai, X., Maskell, K., Johnson, C.A., Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA, 881 pp. ISO, 2005. Soil quality - Determination of pH. International Organization for Standardization 10390, Available at: https://www.sis.se/api/document/preview/ 905787/. (Accessed on July 23, 2018). Kishchuk, B.E., 2004. Soils of the Ecosystem Management Emulating Natural Disturbance (EMEND) experimental area, northwestern Alberta. Nat. Resour. Can., Can. For. Serv., North. For. Cent., Edmonton, Alberta. Inf. Rep. NOR-X-397, 136. Kishchuk, B.E., Quideau, S., Wang, Y., Prescott, C., 2014. Long-term soil response to variable-retention harvesting in the EMEND (ecosystem management emulating natural disturbance) experiment, northwestern Alberta. Can. J. Soil Sci. 94, 263–279. Kishchuk, B.E., Thiffault, E., Lorente, M., Quideau, S., Keddy, T., Sidders, D., 2015. Decadal soil and stand response to fire, harvest, and salvage-logging disturbances in the western boreal mixedwood forest of Alberta. Canada. Can. J. For. Res. 45, 141–152. Kleja, D.B., Svensson, M., Majdi, H., Jansson, P.E., Langvall, O., Bergkvist, B., Johansson, M.B., Weslien, P., Truusb, L., Lindroth, A., Ågren, G.I., 2008. Pools and fluxes of carbon in three Norway spruce ecosystems along a climatic gradient in Sweden. Biogeochemistry 89, 7–25. Laganière, J., Boča, A., Van Miegroet, H., Paré, D., 2017. A tree species effect on soil that is consistent across the species’ range: the case of aspen and soil carbon in North America. Forests 8, 113–128. Laganière, J., Paré, D., Bergeron, Y., Chen, H.Y.H., Brassard, B.W., Cavard, X., 2013. Stability of soil carbon stocks varies with forest composition in the Canadian boreal biome. Ecosystems 16, 852–865. Lenth, Russell V., 2016. Least-squares means: the R Package lsmeans. J. Stat. Softw. 69, 1–33.

8