Climate impacts of retention forestry in a Swedish boreal pine forest

ABSTRACT

The excessive simplification of forest structure associated with clear-cutting can carry risks for biodiversity. Tree retention is an alternative practice that maintains greater structural diversity, but its effects on climate change impacts relative to conventional harvesting are largely unexplored. By integrating field measurements from 12 forest stands with modelling approaches, we investigated the post-harvest effects on radiative forcing of tree retention in a Swedish boreal pine forest. In the near-term, impacts from carbon fluxes and surface albedo were of the same order of magnitude but with opposite sign (warming for carbon and cooling for albedo), with a net warming effect. In the long-term, the net effect turns to cooling, as the forest becomes a strong carbon sink. Retention had a significant effect on climate in the near-term, where differences between the various retention levels are more evident. At increasing tree retention, warming impacts from carbon fluxes tend to decrease, but cooling contributions from surface albedo are less pronounced. Balancing these effects, we find a net climate warming that increases with tree retention.

KEYWORDS:

• Forest management
• climate impacts
• biodiversity
• tree retention
• biodiversity

EDITED BY:

• Karl-Heinz Erb

1. Introduction

Understanding how intensified human interventions affect both climate regulation services provided by terrestrial ecosystems and biodiversity is essential to the deployment of future sustainable use of forests (Bryan et al., 2015; Edwards et al., 2014). Forest ecosystems are serving a variety of purposes, including maintenance of biodiversity and supply of ecosystem services such as timber production, climate regulation, and recreation (Pohjanmies et al., 2017). There is often a trade-off between these different functions. For instance, increased wood outtakes from forests threaten biodiversity and other ecosystem services (Bouget, Lassauce, & Jonsell, 2012; Krause et al., 2017; Kraxner et al., 2013; Michelsen, 2008). In particular, the simplification of forest structure associated with clear-cutting can perturb natural conditions, thereby carrying risks for maintaining adequate life-supporting ecosystem services (Thompson et al., 2011).

1.1. Retention forestry

On a global basis, clear-cut logging is the most common harvesting method (FAO, 2010). It is likely to become more frequent in the future according to the predicted growing demands for energy production from biomass (Slade, Bauen, & Gross, 2014), especially under stringent climate change mitigation objectives (Lauri et al., 2017). At final harvest, retention of trees can mitigate the biodiversity risks associated with clear-cutting (Fedrowitz et al., 2014). Retention forestry is based on partial cutting practices that leave on site single trees, tree groups and dead wood (both standing and lying). Retention increases structural variation in the developing stand, thereby favoring forest-dwelling flora and fauna (Lassauce, Paillet, Jactel, & Bouget, 2011; Radu, 2006). Retained trees and dead wood helps maintaining connectivity in the managed forest landscape and achieving temporal and spatial continuity of key habitat components for both early and late successional species (Bauhus, Puettmann, & Messier, 2009). Retention forestry can also increase public acceptance of forest harvesting and planning of future forest use, with positive effects on water quality and aesthetical values (Clinton, 2011; McDermott, Cashore, & Kanowski, 2010). On the other hand, tree retention can reduce incomes of forest owners owing to lower timber extraction volumes, because it increases costs per unit of wood harvested (Ranius, Ekvall, Jonsson, & Bostedt, 2005; Santaniello et al., 2016). Today, retention forestry is increasingly applied in northern countries like Finland, Sweden, and Canada, with retention levels varying from 1–3% of the harvested volume up to 40% in some Canadian provinces (Gustafsson et al., 2012). In Sweden, certification standards require the retention of green trees (about 10 per ha) and the creation of dead wood in the form of high-cut stumps or girdled trees (about 3 per ha) during harvest (Simonsson, Gustafsson, & Östlund, 2015).

1.2. Forests and climate

Climate change implications of forest management are typically assessed on the basis of absorption or emission of CO₂ and other greenhouse gases (Hudiburg, Law, Wirth, & Luyssaert, 2011; Mitchell, Harmon, & O’Connell, 2012). However, land use and land cover changes also influence the climate by modifying the biophysical properties of the land surface, with potentially large implications for the local and global climate (Betts, Falloon, Goldewijk, & Ramankutty, 2007; Bonan, 2008; Cherubini, Huang, Hu, Toelle, & Stromman, in press). For example, surface albedo, defined as the ratio of reflected to incoming shortwave radiation, largely varies with vegetation cover and local climate. Forest harvest alters the solar reflective property of the surface, as the canopy clearance reduces the surface masking effects of trees, which is particularly relevant in the areas affected by seasonal snow cover (Anderson et al., 2010; Cherubini, Vezhapparambu, Bogren, Astrup, & Strømman, 2017). This tends to increase surface albedo, meaning that a larger fraction of incoming solar radiation is reflected to the atmosphere and not trapped to warm the biosphere (Bonan, 2008). When a forest grows, surface albedo values return to the pre-harvest level.

A change in albedo from afforestation, deforestation, or forest management are found to counteract the effects of the carbon exchanges on the global radiation balance (Arora & Montenegro, 2011; Luyssaert et al., 2014; Mykleby, Snyder, & Twine, 2017; Naudts et al., 2016). In general, climate impacts from carbon fluxes usually dominate at low latitudes, while biophysical contributions from changes in albedo are stronger at high latitudes, where many studies in boreal forest find a nearly climate neutral or cooling effect following clearing (Arora & Montenegro, 2011; Bala et al., 2007; Bathiany, Claussen, Brovkin, Raddatz, & Gayler, 2010; Mykleby et al., 2017). Impacts of forestry operations on climate are therefore to be assessed beyond a simple carbon accounting perspective (Cherubini et al., 2016; Jackson et al., 2008; Mahmood et al., 2014; Zhao & Jackson, 2014).

1.3. Aim and scope

To the best of our knowledge, no studies have explored the influence of retention forestry on global warming. In this work, we integrate impacts from post-harvest carbon and albedo dynamics across a tree retention gradient by coupling field data with modelling approaches. Our hypothesis is that tree retention can have two contrasting effects on climate. On one hand, it reduces the cooling benefits from albedo, because leaving standing vegetation increases snow-masking effects in winter and spring, thereby dampening the post-harvest albedo increase when compared to clear-cut sites. On the other hand, less carbon is harvested and oxidized, thus leading to reduced warming. How the two effects will play out together remains an open question. Hence, the main objectives of our analysis are: i) to quantify the effects of tree retention on climate impacts from carbon fluxes and albedo dynamics, and ii) use this quantification to investigate correlations between trees retained and net climate impacts of harvesting operations.

To address these research questions, we focused on the albedo and carbon dynamics after harvesting 12 pine-dominated forest stands in central Sweden where a tree retention experiment was conducted. We quantified impacts in terms of radiative forcing (RF), a measure of the perturbation to the global radiation budget conventionally interpreted as an approximation of the impacts to average surface temperature. Greenhouse gas (GHG) emissions from harvest operations were also included in the analysis.

2. Methodology

2.1. Study area

The study area, Effaråsen, is located in the province of Dalarna in the southern boreal vegetation zone of Sweden (60°58’29’’N, 14°01’55’’E, mean altitude 385 m). Since 2012, Effaråsen is hosting a retention experiment, in which trees and dead wood, old and newly created, are retained in different amounts across 12 forest stands (Santaniello et al., 2016). The study area has an extension of about 140 ha. Before harvest, it was dominated by 120–140 years old pine trees (Pinus sylvestris), including some older trees. The ground vegetation is dominated by dwarf shrubs (Vaccinium vitis-idaea and V. myrtillus) and lichens (mainly Cladonia spp.). The mean site productivity varies between 2.4 and 4.3 m³ ha⁻¹ year ⁻¹ (mean 3.1 m³ ha⁻¹ year⁻¹). The mean number of stems per ha was 407 (range 212–539), the mean tree height 18 m, and the mean diameter at breast height 26 cm.

2.2. Collection of field data

The retention levels in the harvested stands approximately aimed at about 3% (very low retention), 10% (low retention), 30% (mid retention) and 50% (high retention) of trees and created dead wood. There were 12 stands in total, three per each retention level. The stands varied in terms of forest volumes and number of standing trees before harvest, and retention levels were relative to the pre-harvest situation in each stand. This means that stands within the same retention category can differ in terms of both extraction volumes and number of trees retained. Stands also differed in the site productivity. There were six sites with lower site index (TM3, ETO3, TM10, ETV10, ETV30, and EFF50) and six with high site index (Santaniello et al., 2017). See Table 1 for specific characteristics of the stands.

Table 1. Characteristics of the forest stands, including volumes and carbon content of the forest before and after harvest, productive machine hours (as a sum of harvester and forwarder) and emissions of CO₂, CH₄ and N₂O from harvest operations. Retention category refers to the % of trees left in the stand after harvest. Volumes refer to stem under bark. Total C include carbon stocks in vegetation (below and above ground) and soil.

Stand name	Retention category	Volume before harvestm^3/ha	Total C before harvesttonC/ha	Volume retainedm^3/ha	Trees retainedN/ha	Wood harvestedtonC/ha	Productive machine hourshr/ha	CO2kg/ha	CH4kg/ha	N2Okg/ha
TM-3	Very low	110	77	3	5	21	9.20	736	11.1	28.2
ETO-3	Very low	133	83	6	20	6	10.7	227	3.41	8.71
DS-3	Very low	40	64	4	12	25	2.84	859	12.9	33.0
TM-10	Low	86	74	11	29	15	6.24	499	7.50	19.1
ETV-10	Low	74	72	9	23	12	5.19	415	6.24	15.9
KS-10	Low	142	87	11	32	25	10.8	867	13.0	33.2
ETV-30	Medium	71	72	17	41	9	4.81	384	5.78	14.7
ETO-30	Medium	81	76	19	67	11	6.37	509	7.65	19.5
DS-30	Medium	135	87	18	65	22	12.1	964	14.5	37.0
ETO-50	High	99	83	33	92	11	5.83	485	7.29	18.6
EFF-50	High	97	82	31	110	11	6.07	466	7.01	17.9
KS-50	High	102	91	43	165	11	6.22	498	7.48	19.1

In each stand, the retained trees at harvest were divided into four types, in similar proportions: intact green living trees, girdled trees (i.e. bark stripped along a part of the stem by the harvester machine), high-cut trees (i.e. leaving about 3 m high stumps), and felled trees left on the ground.

We surveyed all living and dying trees and dead wood objects (diameter≥10 cm and length>50 cm) before and after harvest in a 100 × 100 m square plot randomly placed within the center of each stand. We measured tree diameter at breast height (DBH) for every tree and the height of a subset of trees. A function predicting height from the measured DBH was fitted with data and used to obtain DBH and height for every tree. The H100 (dominant height at 100 years of age) site index for each stand was estimated using standard methods based on site characteristics (altitude, latitude, field and bottom layer vegetation, water availability, and soil texture). Further details are available in (Santaniello et al., 2017).

Data on fuel consumption during harvest operations were gathered from onboard computers of the harvester (John Deere 1470, engine power: 160 kW) and forwarder (John Deere 1710, engine power: 190 kW). Consumption was expressed in productive machine hours (G15 minutes per m³ harvested wood without bark, i.e. including machine stops shorter than 15 min.). Consumption data included processing time and driving time within the stands and from the stands to the road. Forwarding distances ranged from 150 m up to 650 m. Time (G15 min) spent by the forest machines to cut 1 m³ of wood under bark was converted first in m³ wood per hectare and then in liters of fuel consumed. We estimated the amount (kg/ha) of carbon dioxide (CO₂), methane (CH₄) and nitrous oxide (N₂O) emitted by fuel combustion using the Ecoinvent database version 3 (Wernet et al., 2016).

Local climate variables are processed as multi-year (from 1948 to present) monthly means (see Supplementary Figure S1). Surface (2 m) air temperature was taken from the historical database of the Swedish Meteorological Institute (SMHI) (Malung station, close to the sampled forest stands), and snow water equivalents (SWE) were gathered from the Nasa Global Land Data Assimilation System Version 2 (GLDAS-2) (Rodell et al., 2004).

2.3. Post-harvest volume and carbon dynamics

We predicted the post-harvest forest developments in the 12 studied stands for 100 years using Heureka, a decision support system for forestry in Sweden (Wikström et al., 2011). Heureka contains deterministic models that can be applied to single stands, forest holdings, or entire landscapes and regions. The post-harvest growth of the trees, including those retained, is simulated for five-year intervals using growth functions described and validated by (Fahlvik, Elfving, & Wikström, 2014). The growth functions as well as functions for tree mortality are based on data from the Swedish National Forest Inventory. Tree growth is species specific and determined by tree, stand and site characteristics. Tree mortality is simulated with functions by (Elfving, 2014) and ingrowth of new saplings in the stands according to (Wikberg, 2004). The carbon stock of dead wood is modeled considering the influx and decay of dead wood. The influx is modeled according to the volume of dead organic matter from tree mortality (Fridman & Ståhl, 2001). The decay of dead organic matter is modeled using an exponential decay rate where the decomposition rate decreases with increased diameter of the dead wood (Ortiz, Hammar, Ahlgren, Hansson, & Stendahl, 2016). Soil carbon development is simulated with the Q-model (Rolff & Ågren, 1999), which simulates the decomposition of soil organic matter based on the continuous quality theory (Agren & Bosatta, 1998). Then, the soil organic matter is assumed to change continuously over time due to decomposition. Carbon in living trees is calculated using biomass expansion factors for components above ground (stem, bark, branches, and needles) (Petersson et al., 2012) and below ground (stumps and roots) (Petersson & Ståhl, 2006). These functions estimate the tree biomass fractions (foliage, branches, stems, stumps, and coarse roots) using information on tree species, stem size, and the allometric proportion with the tree basal area (Eliasson, Svensson, Olsson, & Ågren, 2013).

We created three series of post-harvest forest growth scenarios for each of the 12 stands under different assumptions regarding tree competition for growth and tree density in the stands. Changes in height and diameter of single trees over time are affected by tree competition and modeled using a function from (Elfving & Kiviste, 1997). The area occupied by the retention trees and designated as a competition zone was estimated considering a buffer zone from the central point of single retention trees and trees in small groups. As a sensitivity analysis, we assumed three alternatives to model tree growth competition with radius equal to 0 m (Buffer 0, i.e. no competition zone), 4 m (Buffer 1), and 6 m (Buffer 2) (Öhman, Edenius, & Mikusiński, 2011). The assumed width of the buffer zones was based on the fact that the height and growth of young trees of Scots pine is reduced within a few meters around single retention trees (Elfving & Jakobsson, 2006). We assumed that in the buffer zones, younger trees are constrained by competition while there are no effects beyond the radius. Within the buffer zone, only existing trees are growing, while regeneration of both planted and naturally regenerated seedlings was assumed outside the zone. Areas for retention and production of the plots were calculated using ArcMap 10.3.1.

In the post-harvest forest stands, 2300 containerized seedlings (genetically improved with 10% gain in growth) were planted per hectare after soil scarification. After tree planting at year zero, the trees were assumed to be left growing and dying according to the natural development of a pine dominated stand in central Sweden for 100 years. We used this regrowth scenario, e.g. excluding thinning, because its simplicity allows to better single out the effects of tree retention.

We produced 36 scenarios of post-harvest carbon dynamics, made by three buffer zones per each of the 12 stands. Model outputs were in the form of stem growth in m³ under bark and total carbon stock dynamics in tons of carbon per hectare. The total carbon stock in each forest stand changed over time as a result of the carbon gradually sequestered by the replanted and retained trees and the carbon emitted from oxidation of dead organic matter below and above ground. Simulated dynamics were for 100 years after harvest with a 5-year interval time step, transformed into an annual time step via interpolation. Net ecosystem productivity (NEP), defined as the difference between carbon sequestration in new trees via net primary productivity (NPP) and CO₂ emissions from heterotrophic respiration (Rh), was computed by performing the first order derivative of the post-harvest carbon stock dynamics. The resulting NEP profiles (Supplementary Figure S2) represented the net exchange of CO₂ between the land and the atmosphere that is commonly used in earth surface studies to describe the ecosystem carbon response to forest disturbances (Goulden et al., 2011).

2.4. Post-harvest albedo dynamics

Monthly-mean albedo dynamics were predicted using a simultaneous unmixing and non-linear regression model (Hu, Cherubini, Vezhapparambu, & Strømman, 2017) based on multi-year satellite retrievals of MODIS surface albedo (MCD43A3) combined with high resolution forest inventory databases and meteorological records (Cherubini et al., 2017). The model has specific functions and parameters for different tree species (pine, spruce, and birch), and age or volume can be used as forest structure variables. The predictions from this model were tested against satellite retrievals, and showed reasonably good statistical scores (R² = 0.713, N ≈ 380,000) (Hu et al., 2017). Pre-harvest stem volumes (in m³ per ha) from Table 1 were used to compute the monthly mean pre-harvest albedo values in the forest stands. Post-harvest stem volume dynamics from Heureka were used to predict post-harvest albedo over the timeframe of the analysis. Changes in surface albedo attributed to the harvest event were then computed by the difference between the pre- and post-harvest predicted albedo values (see Supplementary Figure S3, where the computed changes in surface albedo for all the stands of the analysis at 1 and 50 years after harvest are shown).

2.5. Radiative forcing (RF)

The temporal change in atmospheric CO₂ concentration f(t) from the specific NEP profiles and CO₂ emissions from the harvested wood (h_s, assumed to be instantaneously oxidized) was computed for each stand s and buffer zone b after the integration with the global carbon cycle through a mathematical convolution:

(1)

Where y(t) is the impulse response function to a CO₂ emission pulse, simulated using a multi-model mean (Joos et al., 2013) in line with the last IPCC assessment report (Myhre et al., 2013). This value was then translated into RF using the radiative efficiency of CO₂, which is the marginal increase in RF following a unit (kg) increase in the atmospheric abundance of CO₂ (relative to the average atmospheric CO₂ concentration in 2010 of 389 ppmv) (Myhre et al., 2013). The instantaneous RF was computed under the common assumption that for sufficiently small emissions and approximately constant background conditions the radiative efficiency can be approximated as time-invariant (Myhre, Highwood, Shine, & Stordal, 1998).

Impacts from surface albedo changes are frequently compared with those of greenhouse gases in terms of radiative forcing (Betts, 2000; Betts et al., 2007; Lutz & Howarth, 2014; Randerson et al., 2006). We used a simplified 1-layer atmospheric transfer model where a time-varying local change in monthly-mean surface albedo (Δα(m,t)) is translated to a change in global radiative forcing at the top of the atmosphere (RF^α) after combination with the monthly-mean incoming solar radiation at the surface (R(m)) and a global average atmospheric transmittance parameter (T) (Bright et al., 2014; Cherubini, Bright, & Strømman, 2012):

(2)

where A_Earth is the area of the Earth. An average of the values of RF per month m provided the annual mean radiative forcing for the period of the analysis per each forest stand and buffer zone simulation. Spatially explicit mean incoming solar radiation R was gathered per month and latitudinal/longitudinal degree as the 22-year mean of the downwelling solar radiation flux at surface from the NASA Power Project. A multi-year mean was here used to reduce the noise from inter-annual variability. The atmospheric transmittance parameter was assumed to be constant (0.874) (Lenton & Vaughan, 2009).

3. Results

3.1. Radiative forcing from the retention levels

In all stands, the change in surface albedo associated with the harvest event caused a short-term radiative forcing of the same order of magnitude as the forcing from the corresponding CO₂ fluxes (Figure 1). Dynamics of radiative forcing impacts were highly time dependent. For the first couple of decades, impacts from carbon fluxes were positive (i.e. warming) and from surface albedo negative (i.e. cooling). Thereafter, the sign switched and contributions from carbon largely overwhelmed those from surface albedo. Carbon impacts turned to negative because the forest gradually becomes a strong carbon sink as it cumulatively sequesters more carbon than the one present before harvest due to the improved tree growth. This means that forest stands can achieve higher standing volumes per hectare than pre-harvest levels, which on the other hand makes albedo turns to positive values after a few years (as the forest is denser).

Figure 1. Radiative forcing following forest harvest across different tree retention levels. Results include contributions to radiative forcing from carbon fluxes and albedo dynamics and aggregate the outcomes of the individual 12 forest stands of the experiment into retention category very low (a), low (b), medium (mid) (c) and high (d). For each forest stand, three buffer zone simulations are considered for modelling post-harvest dynamics, representative of tree growth competition with radius equal to 0 m (Buffer 0), 4 m (Buffer 1), and 6 m (Buffer 2). The average of these simulations in the same retention category are shown with solid lines, and the coloured areas indicate the minimum and maximum values.

In Figure 1, results are aggregated per retention category, so that the shaded areas represent the range given by the stands within the specific retention level and the respective simulations for different buffer zones. Post-harvest growth rates were larger in the simulations with lower tree growth competition (Buffer 0), and gradually decreased at increasing growth competition (from Buffer 1 to Buffer 2) (see Supplementary Figure S2). However, differences among buffer zones was smaller if compared to the overall carbon dynamics, because the effects of growth competition were overwhelmed by the improved growth rates of the regenerating trees.

In absolute terms, impacts from carbon tended to decrease at increasing tree retention, owing to lower fraction of carbon emitted to the atmosphere and dead organic matter created from the harvest event. Radiative forcing impacts start at about 1 × 10⁻¹⁰ W m⁻² ha⁻¹ in the case of low and very low tree retention, and at about 0.75 × 10⁻¹⁰ W m⁻² ha⁻¹ for high retention (Figure 1). On the other hand, impacts from albedo changes tended to decrease with increasing tree retention because lower fraction of the canopy is cleared and thus average values of albedo change are smaller. They start from about −1.2 × 10⁻¹⁰ W m⁻² ha⁻¹ in the very low retention category and from about −0.5 × 10⁻¹⁰ W m⁻² ha⁻¹ for high retention.

The net effects on radiative forcing of carbon and albedo for the different retention categories are shown in Figure 2. In the short-term, average net radiative forcing impacts follow a bell-like shape curve and are positive, except for the very low retention case where impacts are negative for the first few years, and increase with tree retention (Figure 2(a)). The peak in net radiative forcing is about 2.3, 1.5, 1.2 and 0.8 W m⁻² per ha for the high, medium, low, and very low retention level, respectively, with the timing of the peak that gradually occurs later at reduced tree retention. This means that high tree retention causes the highest positive peak (warming) in radiative forcing. In general, tree retention increased radiative forcing relative to the case with very low retention, except for a time interval between 20 and 50 years after harvest where the case with mid retention shows the lowest impacts (Figure 2(b)). In the long-term, differences between retention levels depend on the growth rate and cumulative carbon sequestration in the new trees, which is larger in the stands with very low retention. As shown in Figure 2(b), the very low retention case had the smallest average radiative forcing in the long term. In these stands, planting of new trees achieves higher carbon sequestration levels than the previous forest successions. In the cases with higher retention, older trees that are close to saturation are retained and they make less favorable growing conditions for the new trees.

Figure 2. Net radiative forcing results after combinations of contributions from CO₂ fluxes and surface albedo changes. (a) comparison of net impacts for the average results in the four retention categories; (b) comparison of net radiative forcing in relative terms, using the results of the very low retention category as a benchmark.

3.2. Correlation between RF and tree retention

We find a decreasing trend in radiative forcing per hectare from both carbon and albedo at increasing number of trees and fraction of m³ retained (Figure 3). Ten years after harvest, carbon radiative forcing ranges from approximately 1.5 × 10⁻¹⁰ W m⁻² ha⁻¹ when about 5 trees were retained (and < 5% of volume retained) to about one third when about 165 trees were retained (and > 40–45% of volume retained). Impacts from albedo changes range from about −1.25 x 10⁻¹⁰ W m⁻² ha⁻¹ to about −0.3 x 10⁻¹⁰ W m⁻² ha⁻¹ at increasing number of trees and volume retained. Because impacts from carbon are larger than those from albedo, the net results for all the stands typically lie on the positive side of the domain. Impacts 100 years after harvest were largely insensitive to retention levels, because over the long-term the influence of different levels of tree retention is overwhelmed by the improved growth rates of the replanted trees.

Figure 3. Relationships between radiative forcing per hectare at 10 (a) and 100 years (b) after harvest and key characteristics of the forest stands, such as number of trees retained and fraction of m³ retained. Results (Buffer 1 only) show the net effects and the single contributions from carbon and albedo.

The radiative forcing per hectare normalized to the unit of wood harvested is shown in Figure 4. Stands with higher retention levels had generally higher normalized impacts than those with lower retention (Figure 4). In a couple of cases net results are negative (DS-3 and ETV-10), and impacts can increase more than six times going from a very low retention level (about 0.4 × 10⁻¹² W m⁻² ha⁻¹ tC⁻¹ for TM3) to a high retention level (about 2.6 × 10⁻¹² W m⁻² ha⁻¹ tC⁻¹ for KS-50). The trend in net results has some fluctuations across the retention gradient, which can be explained by the heterogeneities of the stands. The trend is clearer for stands with medium or high retention levels, which had more homogeneous pre-harvest characteristics. When the normalized results are ranked as a function of the number of trees retained (Figure 4(b)), the increasing trend of climate impacts with tree retention becomes clearer: the more trees are retained at harvest the higher is the warming impact attributed to the harvested wood per ha.

Figure 4. Radiative forcing per ha normalized to the unit of wood harvested (expressed as tons of carbon, tC) 10 years after harvest (Buffer 1 only). (a) normalized contributions from carbon and albedo and net results for the individual forest stands; (b) normalized net radiative forcing as a function of number of trees retained.

3.3. GHGs from harvest operations

Radiative forcing impacts from emissions of GHGs, mainly CO₂, with smaller contributions (≤ 1%) from CH₄ and N₂O, from harvest operations are about one order of magnitude smaller than those from carbon and albedo (Figure 5). Low retention tends to have lower fuel consumption and hence lower GHG emissions than high retention. Harvest intensities usually increase with the amount of carbon that was harvested rather than with the number of tree retained (Figure 5(b)). In general, we did not observe a strong correlation between trees retained and radiative forcing from harvest emissions. Harvest operations are dependent on other factors than tree retention, such as morphology of the terrain, distribution of the standing trees, and distance to forest road, which vary for the different stands.

Figure 5. Radiative forcing impacts for each of the individual forest stands of GHG emissions from harvest operations (a) and as a function of the corresponding amount of wood harvested (b).

4. Discussion

Balancing all the effects together, we find that tree retention generally increases net warming impacts per ha of land (Figure 3(a)) and harvested wood (Figure 4(b)). However, radiative forcing perturbations from forest harvest with different levels of tree retentions have complex and heterogeneous dynamics. Retaining trees reduces impacts from carbon, but tree retention also reduces cooling contributions from surface albedo changes, as a closer canopy is retained relative to clear-cut. Referring back to our initial hypothesis, we see that the net effect of these two mechanisms changes over time. In the first years we find an average warming contribution, as the cooling effects from albedo are overwhelmed by the warming contributions from carbon. However, a couple of decades after harvest, the situation changes and the forcing switches sign. Radiative forcing from CO₂ turns to deep negative values owing to the more favorable growing conditions of the forest succession, as larger amounts of carbon accumulates in the forest than those stored before harvest. Impacts from albedo become positive due to higher tree density and stabilize after about 60 years. The net effect thus turns to cooling as it is dominated by the large carbon sequestration fluxes.

4.1. Interpretation of RF from carbon dynamics

The near-term increase in radiative forcing from carbon is due to the time lag between CO₂ emission and removal fluxes that caused a short-term increase in atmospheric CO₂ concentration. This is a trend observed in previous studies of carbon dynamics in forest successions (Bernier & Paré, 2013; Cherubini, Gasser, Bright, Ciais, & Stromman, 2014; O’Halloran et al., 2012). Soon after harvest, carbon impacts are high due to both emissions from harvested wood, whose carbon is assumed to be instantaneously released to the atmosphere, and decay of dead organic matter from the forest. Forest stands are a source of carbon for some years after harvest, especially if they have low retention levels, because CO₂ emissions from heterotrophic respiration exceed carbon sequestration in trees via NPP, and NEP is therefore negative. Once residues have decomposed and NPP increases, NEP becomes positive, and the forest ecosystem acts as a net carbon sink (Goulden et al., 2011; Harmon, Bond-Lamberty, Tang, & Vargas, 2011). In our study, the transition from carbon source to carbon sink is observed within a time period up to 15 years after harvest. This timing of the NEP profiles is consistent with other studies that measured post-harvest carbon dynamics in boreal forest via chronosequences based on flux tower data (Amiro et al., 2010; Zha et al., 2009). These dynamics are evident in Figure 1. After the initial period, RF from carbon declined as the sink of the growing trees progressively sequesters carbon from the atmosphere, and became negative because the stand sequestered more total carbon than what was standing before harvest. This makes the forest management activities carbon negative from the medium to the long term, as higher levels of carbon are stored per unit of land. The decrease in the breadth of the shaded area at increasing retention category is mainly due to the homogeneity of the forest stands with high tree retention. As mentioned above, stands with very low retention have more differences in the pre-harvest standing volumes, which makes their results fall in a wider range.

4.2. Interpretation of RF from albedo dynamics

Compared to carbon, impacts from albedo were of opposite sign and had the same order of magnitude in the short-term, but they were significantly smaller in the long-term. The removal of the canopy increased surface albedo creating a cooling effect that decreased only slightly for the first couple of decades. Forcing from albedo then decayed with the gradual regrowth of the trees and the associated progressive closure of the canopy. Albedo differences were larger in winter than in summer. Supplementary Figure S3 shows that annual means of albedo changes soon after harvest were between 0.04 and 0.15 in January and 0.04 and 0.01 in July.

Albedo changes were larger in the cases of lower tree retention. A special case was DS-3, where there was little standing volume before harvest (40 m³), and the associated change in albedo after harvest was thus smaller than that of the other stands in the same retention category (Figure S3a). Similarly, initial albedo changes were larger for the cases in which standing volume before harvest was high, and tree retention was low, like ETO-3 (133 m³), TM3 (110 m³), KS-10 (142 m³). On the other hand, they were small when the stands had low standing volume (as DS-3) or high tree retention (e.g. KS-50 and ETO-50).

In all cases, albedo changes decreased and turned to positive after about 20–50 years (Supplementary Figure S4), when the average albedo values were smaller than those at pre-harvest levels (see Supplementary Figure S3b for monthly mean albedo values in year 50 after harvest). DS-3 resulted in the largest difference in albedo in the long term. This is a highly productive site that achieved a much higher level of standing volume in relation to pre-harvest values.

Radiative forcing from albedo reached the maximum cooling in March, where there was the best combination of land surface albedo values (as there is usually still snow on the ground) and incoming solar radiation (see Supplementary Figure S5 for the monthly mean radiative forcing from surface albedo change the first year after harvest and after 100 years).

4.3. Uncertainties and limitations

The study has limitations pertaining to the modelling approaches used and to the field measurements undertaken to produce the data used in our analysis.

Forest characteristics are heterogeneous across the stands. The retention categories were based on the relative number of trees left per stand, but stands in the same retention category had varying numbers of trees and carbon stocks, which make the volumes of harvested wood and number of trees retained differ (Table 1). This was especially true for the stands in the very low retention category. For example, TM-3 and DS-3 had large differences in volume before harvest, 110 against 40 m³/ha, and apparently lower differences in pre-harvest total carbon stock. The average soil carbon content of these stands was about 56 tC/ha, and 5 m³ of wood corresponds to about 1 ton of carbon. In the stand TM-3, 110 m³ thus corresponded to about 21 tC/ha, giving a total carbon content (soil + vegetation) before harvest of 77 tC/ha. In the stand DS-3, the total carbon content before harvest was 64 tC/ha (56 tC/ha from the soil and 8 tC/ha from the vegetation). Similarly, 20 retained trees only made 6 m³ in ETO3 because it had a mean volume of about 0.3 m³, whereas the 5 retained trees in TM3 had larger average volume (about 0.6) and resulted in 3 m³/ha. In general, DS-3 represented a special case because the stand had relatively little and sparse standing volume before harvest. Homogeneity of the stands generally increased for those within the medium and high retention category.

The modelling of post-harvest carbon dynamics and tree growth is based on empirical data from more than 14 000 National Forest Inventory (NFI) plots. The starting values for soil carbon in the post-harvest simulations were not empirical data from the experiment plots. They were based on the Q-model that relies on empirical data from the Swedish Forest Soil Inventory and uses input data about site factors (altitude, latitude, soil texture, humidity and vegetation cover) and stand data from the experiment plots. In this study, we focused on the changes rather than the size of the soil carbon pool, thereby minimizing the uncertainty associated with estimates of soil carbon contents.

Retained trees can also show divergent patterns in mortality and growth rates. Mortality is found to increase in cases of wet soils and high wind exposure, and decreases at increasing tree volume in retention groups (Hallinger, Johansson, Schmalholz, Sjöberg, & Ranius, 2016), whereas volume growth response to partial cutting can increase (Bose, Brais, & Harvey, 2014). These effects are inherently difficult to model and their influence is expected to gradually increase under high tree retention levels. However, their impact on the results is dampened by the opposite direction that these effects have on the carbon fluxes.

The model used to simulate pre- and post-harvest albedo dynamics was produced from an extensive dataset of albedo satellite retrievals integrated with high resolution land cover maps and meteorological records. Originally derived for a forest area at the border between Norway and Sweden, its transportability tests did not report significant biases when compared against observations (Hu et al., 2017).

Since the radiation budget is the fundamental driver of the climate system, radiative forcing is a useful indicator of the relative importance of the different mechanisms on global temperature changes, and it is frequently used to assess different forestland perturbations on climate (Andrews, Betts, Booth, Jones, & Jones, 2016; Betts et al., 2007; Lohila et al., 2010; O’Halloran et al., 2012; Randerson et al., 2006). To compute radiative forcing from surface albedo changes, we applied a simplified radiative forcing model used in previous studies (Bright et al., 2014; Caiazzo et al., 2014; Cherubini et al., 2012), which perform reasonably well when compared with more sophisticated models (Bright & Kvalevåg, 2013). The approach used to compute radiative forcing from carbon fluxes aligns with the IPCC standard protocol for quantification of metrics and climate impacts (Myhre et al., 2013), where forcing impacts are computed under a constant background climate, and possible feedbacks of climate change into the terrestrial carbon cycle are excluded.

Future climate change and increasing atmospheric CO₂ concentrations are expected to cause substantial changes in vegetation structure and functions. Understanding the net effect of climate feedbacks on terrestrial vegetation is complex because of the variety of factors involved. Changes in NPP, Rh and NEP occur primarily through enhanced growth because of changes in climate conditions and CO₂ fertilization effects, nitrogen availability, faster decay rates of residues, and higher soil respiration (Gerber, Hedin, Keel, Pacala, & Shevliakova, 2013; Loudermilk et al., 2013). However, the overall effect of projected future climate change on terrestrial vegetation has large heterogeneity in terms of the magnitude and sign of the net change in NEP, and the response of the terrestrial carbon cycle to climate change is one of the largest sources of uncertainty affecting future climate change projections (Friend et al., 2014).

A changing climate will also influence radiative forcing contributions from surface albedo (Brovkin et al., 2013), owing to the predicted average temperature rise that can reduce snowfalls and thereby surface albedo. The possible influence of this feedback on the conclusions of this study is reduced because the major forcing from albedo is exerted in the first couple of decades after harvest, and the albedo feedback is to some extent compensated by the decrease in the marginal radiative efficiency of CO₂ (Bright et al., 2014).

The concept of radiative forcing has some limitations itself. For example, climate sensitivity to a given level of forcing can vary depending on the characteristics of the forcing and the latitude at which it operates (Hansen et al., 2005; Shindell, Faluvegi, Rotstayn, & Milly, 2015), and when translating into climate impacts additive RF from land cover changes and GHGs there can be differences in temperature and precipitation responses (Jones, Collins, & Torn, 2013). Radiative forcing cannot be used to quantify all mechanisms of climatic perturbation from land use, such as those which act directly via surface moisture fluxes (Davin & de Noblet-Ducoudré, 2010; Pielke et al., 2011). The latter directly influence near-surface temperature and local climate, and, along with impacts from trace emissions of other climate-active compounds, affect the climate system though complex mechanisms that are highly site specific and difficult to quantify (Unger, 2014; West, Narisma, Barford, Kucharik, & Foley, 2011). However, changes in surface albedo are found to dominate other biophysical effects at high latitudes, where carbon and albedo are the major drivers for climate change (Bala et al., 2007; Claussen, Brovkin, & Ganopolski, 2001; Davin & de Noblet-Ducoudré, 2010).

5. Conclusions

This study is the first to assess impacts from tree retention on climate warming from forest harvesting. The analysis shows the importance to include contributions from surface albedo changes, as they are of the same order of magnitude as impacts from carbon fluxes in the short-term. Tree retention affects radiative forcing impacts from both carbon and albedo dynamics. We found that increasing the number of retained trees in managed forests increases climate warming.

Future studies can deepen on the resolution and better target the individuation of a sensitivity threshold for the retention level, specifically based on the interactions between retained trees required to maintain species richness and climate. The identification of an optimum threshold for the number of trees to be left standing, or type and quantity of dead wood to be retained, would allow the elaboration of case-specific forest management strategies that can co-deliver in terms of the dual target of climate change mitigation and biodiversity conservation.

Overall, quantifying climate change impacts from forest harvest disturbances is complex, as the outcomes vary greatly across space, time and disturbance intensity. Improving our understanding of the processes regulating forest ecosystems functioning and the climate is key to design the best sustainable forest management strategies under the multiple roles that forests can play. Existing models and approaches can be tweaked and integrated to achieve interdisciplinary research involving forest ecologists, climate and environmental scientists to investigate trade-offs and synergies in terms of resource supply, climate change mitigation and biodiversity conservation.

Disclosure statement

No potential conflict of interest was reported by the authors.

Supplementary material

Supplemental data for this article can be accessed here.

Additional information

Funding

The work of F.C., X.H., and A.H.S. was supported by the Norges Forskningsråd [244074].