Net soil carbon balance in afforested peatlands and separating autotrophic and heterotrophic soil CO₂ effluxes

Peatlands are a significant global carbon (C) store, which can be compromised by drainage and afforestation. Quantifying the rate of C loss from peat soils under forestry is challenging, as soil CO₂ efflux includes both CO₂ produced from heterotrophic peat decomposition and CO₂ produced by tree roots and associated fungal networks (autotrophic respiration). We experimentally terminated autotrophic below-ground respiration in replicated forest plots by cutting through all living tree roots (trenching) and measured soil surface CO₂ flux, litter input, litter decay rate, and soil temperature and moisture over 2 years. Decomposition of cut roots was measured and CO₂ fluxes were corrected for this, which resulted in a large change in the fraction heterotrophic : autotrophic flux, suggesting that even 2 years after trenching decaying root biomass makes significant contributions to the CO₂ flux. Annual peat decomposition (heterotrophic CO₂ flux) was 115 ± 16 g C m⁻² yr⁻¹, representing ca. 40 % of total soil respiration. Decomposition of needle litter is accelerated in the presence of an active rhizosphere, indicating a priming effect by labile C inputs from roots. This suggests that our estimates of peat mineralization in our trenched plots are conservative and underestimate overall rates of peat C loss. Considering also input of litter from trees, our results indicate that the soils in these 30-year-old drained and afforested peatlands are a net sink for C, since substantially more C enters the soil organic matter than is decomposed heterotrophically. This study does not account for fluvial C fluxes, which represent a small flux compared to the CO₂ soil efflux; further, root litter and exudate deposition could be a significant C source that is only partially sampled by our approach, adding to these plantations being a potential carbon sink. However, the C balance for these soils should be taken over the lifespan of the trees, in order to determine if the soils under these drained and afforested peatlands are a sustained sink of C or become a net source over longer periods of forestry.

1 Introduction

Large peatland areas in the boreal and temperate zone have been drained and afforested in western Europe, especially in the UK, Ireland and the Fennoscandia region, with conifers replacing native peatland vegetation (Andersen et al., 2016). In the UK alone, > 800 000 ha (approximately 20 %) was affected by this land-use change (Anderson et al., 2016). As well as causing habitat loss, drainage and afforestation of peatlands influence peatland hydrology, biogeochemistry and carbon (C) storage. In Fennoscandia, drainage to improve tree growth on naturally treed peatlands has been shown to cause significant changes in decomposition and greenhouse gas (GHG) fluxes (Ojanen et al., 2013). However, the effects of drainage and afforestation on temperate peatlands that were previously dominated by less-statured vegetation (e.g. blanket bogs, moorland, heathland) are more uncertain (Sloan et al., 2018). While changes in hydrology and soil redox potential are anticipated to accelerate soil C loss and alter the composition of GHG emissions, there are very few quantitative data on how this shift in land use changes soil C dynamics and GHG emissions in temperate peatlands (Hermans et al., 2018). Since Northern Hemisphere peatlands are estimated to store a third of global terrestrial carbon (Gorham, 2010), changes linked to drainage and afforestation of temperate peatlands could have significant impacts on regional C dynamics and GHG balances.

There is a broad consensus that peatland drainage accelerates the loss of endogenous peat C stocks, but the impacts of drainage and afforestation on total soil C stocks and soil respiration are less certain (Hargreaves et al., 2003; Mayer et al., 2020; Simola et al., 2012; Vanguelova et al., 2019; Zerva and Mencuccini, 2005). Afforestation of previously open peatlands in the British Isles (e.g. Mazzola et al., 2021) and Fennoscandia (e.g. Tolvanen et al., 2020) differs from other types of drained peatlands because this land-use change involves a wholesale shift in the functional composition of the plant community (i.e. replacement of short-statured vegetation with trees), leading to potential interactions or synergistic effects arising from changes to both soil hydrology and plant community composition.

The comparatively high rates of net primary production and larger stature of the trees on drained and afforested peatlands can represent a significant net ecosystem C sink and consequently represent a large source of detrital material to soil (Yamulki et al., 2013). Thus, it is possible that this input of C from more highly productive trees could partially or wholly offset or outstrip some of the losses of C derived from mineralization of the peat. Moreover, changes to soil organic matter (SOM) chemistry due to a shift towards inputs of more recalcitrant and nutrient-poor coniferous litter could further enhance soil C storage in afforested peatlands in the British Isles and Fennoscandia. This is because the trees planted in these regions (e.g. Sitka spruce) produce tissues that are less readily decomposed than the organic matter (OM) generated by the short-statured vegetation that they have replaced (i.e. larger proportion of woody debris generated with higher C : N ratios and greater proportion of recalcitrant compounds like lignin) (Hermans et al., 2018). This may affect soil C stocks by inhibiting or slowing overall rates of decay due to the tree litter's recalcitrance (Liski et al., 2002) but has not been studied in tree plantations on deep peat so far.

In order to close these knowledge gaps, we need to determine the C balance of drained and afforested peatland soils, tracking C inputs from tissue turnover (e.g. litterfall and rhizodeposition), and C losses from SOM decay. However, separating the direct and indirect effects of trees on peat mineralization is methodologically challenging (Subke et al., 2006). Soil CO₂ efflux includes both CO₂ released from peat mineralization (heterotrophic CO₂ flux) and CO₂ produced from root—rhizosphere (i.e. autotrophic) respiration; experimental manipulation of autotrophic C supply to the rhizosphere allows a separation of these two main component fluxes but introduces a number of potential artefacts (see Subke et al., 2006, for a methods overview including method-related assumptions and artefacts). Root trenching, where roots are severed throughout the depth of the rooting zone to create areas of forest floor with no recent C input from trees to roots or associated mycorrhizas, has been used across many forest sites and provides useful insights into relative contributions of autotrophic and heterotrophic soil CO₂ efflux and the respective temporal dynamics. This disruption of root exudation and continuous root turnover could influence peat decay through processes such as microbial priming effects (Kuzyakov, 2006; Subke et al., 2006), where microbial activity (linked with an active rhizosphere) in the soil is stimulated by the addition of easily accessible C from roots (Kuzyakov, 2006). Further experiments are required to evaluate the effects of these other plant-facilitated processes on peat decay to better constrain estimates of heterotrophic activity in soil.

In this study, we experimentally determined the C budget of a drained and afforested peat soil. Carbon inputs were tracked by quantifying litterfall, while C outputs were determined using a trenching technique to partition root–rhizosphere and heterotrophic respiration. In addition, we conducted a root decomposition study in rhizosphere and living root-free (i.e. trenched) soil to constrain artefacts associated with the trenching and in order to determine if rhizosphere-linked processes (e.g. priming effects) facilitated OM decay. We hypothesized that (1) the soils under these forestry plantations are a net C source, (2) soil CO₂ efflux is dominated by autotrophic respiration, and (3) interactions between C supply to the rhizosphere by trees and surface litter decomposition result in higher decomposition rates of the litter due to rhizosphere priming effects.

2 Methods

2.1 Study site

The research took place in RSPB's Forsinard Flows National Nature Reserve in the north of Scotland (58^∘22^′ N, 3^∘53^′ W). Four paired plots were established in the beginning of June 2014 in three separate forestry plantation blocks of identical age containing a mixture, 3 : 5 (on average), of Sitka spruce (Picea sitchensis) and lodgepole pine (Pinus contorta). The plantations were drained, around 30 years old and very dense (about 5000 trees per hectare), with no vascular plant understorey but sporadic patches of moss, predominantly feather moss, e.g. Hypnum jutlandicum, Hylocomium splendens and, in some instances, Sphagnum fallax and S. capillifolium in furrows. The average diameter at breast height (DBH) for Sitka spruce was 13.3 cm and for lodgepole pine 17.9 cm. The average canopy cover was 76.3 % (Smith and Hancock, 2016). Peat depths in these three forestry blocks varied between 30 and 260 cm, with depths at research plots between 137 and 204 cm (Smith et al., 2014). The average annual precipitation in the research area between 1981–2010 was 970.5 mm, with an average maximum air temperature of 11.4 ^∘C and average minimum air temperature of 3.3 ^∘C measured at the Kinbrace weather station approximately 20 km from the plots (location: 58^∘13^′89^′′ N, 3^∘55^′1.2^′′ W; altitude: 103 m above mean sea level; Met Office, 2018). Seasonal averaged water table relative to ground surface is −350 (±20) mm in spring (March–May), −457 (±34) mm in summer (June–August), −404 (±49) mm in autumn (September–November) and −244 (±14) mm in winter (December–February).

2.2 Experimental set-up

Candidate locations for the four paired trenched and control areas in each plot were initially identified at random and soil surface respiration measured. Based on respiration results, two closely matched plots were selected and randomly allocated a treatment (trenching or control). Paired plots were no more than 10 m apart from each other.

The double-mouldboard ploughing of the peat at the time of afforestation created a regular microtopography with low-lying furrows (ca. 1.5 m wide) flanked by high ridges (plough throws; ca. 0.75 m wide) on either side. In between two plough throws, there is an area of ca. 0.50 m width of the original (unploughed) surface (Fig. 1). The height from the bottom of the furrow to the top of the plough throw is on average 0.5 m, and from the original surface to the plough throw it is about 0.15–0.2 m. In general, conifer seedlings were planted on the plough throws because of the improved drainage compared to the original surface. For our study, each plot encompassed the three types of microtopography.

Figure 1 Microtopography in forest plantations, with location of measurements for flux chambers (not to scale).

2.3 Trenching

A 40 cm deep and 30 cm wide trench was cut to just below the main rooting depth of the trees, cutting through all roots present. The trench was double lined using polypropylene gardening cloth and refilled with peat soil in between the two layers of cloth to prevent ingrowth of roots (Fig. 2). The dimensions of each trench plot were about 3.5×1.5 m and included all three microtopographic forms. These dimensions maximized the space between trees but did not represent ground within 60 cm of trees.

Figure 2 (a) Trenching plot dug and all living tree roots are sawn trough. (b) Finished trench plot; trenches double lined with polypropylene gardening cloth and filled in with peat.

2.4 CO₂ measurements

Three pairs of PVC collars of 10 cm height with a diameter of 20 cm were installed to a depth of 3 cm within the three microforms (Fig. 1) of both trenched and control plots. CO₂ measurements were taken using custom-built dark dynamic closed chamber with a height of 5 cm and a diameter of 20 cm, which were placed on the permanent collars for 3 min. Elasticated rubber material placed over the joint of chambers and collars for the duration of the measurement provided an airtight seal. The chamber was connected to an EGM 4 infrared gas analyser (PP-Systems, Amesbury, MA, USA), recording CO₂ concentrations every 4–5 s. Fluxes were calculated from increases in CO₂ concentration within the chamber over 3 min. Measurements were carried out 10 times between August 2014 and July 2016.

2.5 Litter

Six litter traps (0.07 m² each) were located close to each plot, and litter (falling needles and twigs) was collected at each sampling visit. The dry weight of all litter was recorded as area-based averages for each plot separately.

Litter was allowed to fall onto the soil surface within collars for the duration of the experiment. To be able to distinguish peat respiration from litter respiration, surface litter was removed manually from one (always the same) of the two paired collars in each microform before measuring respiration. The litter present in the collar with litter was weighed after a respiration measurement and then placed back in the collar. Litter adjacent to the collar was also collected and weighed in the field, then taken back to the lab, and dried and weighed again to establish the wet-to-dry-mass ratio of litter and calculate litter dry mass within each collar.

2.6 Roots

Root biomass was determined from soil cores (0–20 cm deep and 6.5 cm diameter) taken from each microform in all plots, at the start (June 2014) and end (July 2016) of the experiment. Both dead and living roots from each core were carefully separated and sorted into three root diameter classes: smaller than 2, 2 to 5 and greater than 5 mm. All roots and the root-free soil were dried at 50 ^∘C for 7 d and weighed to establish percentage of roots per gram of soil. Root depth was found to be 25 cm when digging trenches, and root density for 0–20 cm was assumed to be representative for 0–25 cm.

To estimate root decomposition directly, roots were taken from soil collected in each plot, dried (50 ^∘C for 7 d) and separated into the same size classes as described previously. Between 0.36 and 0.69 g of dried root material (separately for each size class) was placed in polyester mesh bags (10×10 cm; mesh size of 0.5 mm) for field incubations. Bags were soaked in water for 2 d prior to field placement to mimic field conditions. Four replicate bags of each size class were buried at 5–10 cm depth in all three microforms in all plots 4 weeks after trenching. To account for any weight loss that may have occurred prior to field incubation, five bags of each size class were taken into the field and not buried but taken back to the lab; the proportional mass loss of litter in these bags was used to correct the initial root mass of all other bags.

One bag per root class per microform was collected from all sites in November 2014, March 2015, July 2015 (except root class >5 mm, since there was not enough material for four bags) and July 2016. After retrieval, bags were dried for 7 d at 50 ^∘C, and root dry mass was recorded.

Root decay was fitted to an exponential decay function:

$$\begin{array}{}\text{(1)}& M_t=M_{\mathrm{0}}{e}^{-kt},\end{array}$$

with M_t the remaining amount of root biomass after collection from the field, M₀ the initial root biomass, and t time and k the decay constant. Data fits were performed separately for root size and microform.

2.7 Soil moisture and temperature

Between June 2014 and July 2016, soil moisture and soil temperature at 5 and 20 cm soil depth were recorded at 30 min intervals in all three microforms in a nearby plot, using 12-bit smart temperature sensors, S-TMB-M002 (Onset Computer Corporation, Bourne, MA, USA), and 10HS soil moisture smart sensors, S-SMD-M005 (Decagon Devices, Inc., Pullman, WA, USA), combined with Onset's smart sensor technology, measuring m³ m⁻³ (volumetric water content) connected to a Hobo micro-station (Onset Computer Corporation, Bourne, MA, USA).

In addition to this, soil temperature (10 cm thermistor) and moisture (HH2 moisture meter, Delta-T Devices, Cambridge, measuring m³ m⁻³ (volumetric water content)) were measured at 5 cm depth next to each collar during sampling. Air temperature was also measured 1 m above the ground during sampling.

2.8 Statistical analyses and flux calculations

Data were analysed using R (R Core Team, 2016) in RStudio version 1.0.136. All CO₂ data were log transformed to meet the criteria of normality, and a linear mixed-effect model (LMM) using the lme4 package in R (Bates et al., 2015) was used to predict CO₂ fluxes based on environmental drivers. All numerical predictors were standardized to 1 standard deviation prior to statistical analyses to allow relative effect sizes of predictors to be compared directly (Nakagawa and Schielzeth, 2010). Model selection was done based on information theory (Burnham and Anderson, 2002). First a model of maximum complexity was built, with interactions between soil moisture, soil temperature, trench treatment, and microform plus interactions between trench treatment, microform and litter treatment, with plot as a random effect (lmer(CO₂ flux ∼ (soil moisture ⋅ soil temperature ⋅ treatment ⋅ microform + treatment ⋅ microform ⋅ litter treatment) + (1|plot))). Then, all possible combinations of this model were identified using the “dredge” function in the MuMIn package (Barton, 2017). Goodness of model fit was assessed with the small-sample-size-corrected Akaike information criterion (AICc), which is calculated using the number of parameters and either the maximum likelihood estimate for the model or the residual sum of squares. “Likelihood” here is a measure of the extent to which a sample provides support for particular values of a parameter in a parametric model. AICc values of different models can be compared, and the model with the lowest AICc is selected as the best approximating model, hereafter called “top model” (Burnham and Anderson, 2002). Any of the models with a delta AICc of less than 2 are considered to be as good as the best model (Richards, 2005). The dredge function also gives a weight to the models it produces, ranging between 0 and 1, with for example a weight of 0.7 meaning that there is a 70 % chance that that model is the best approximating model of the models considered. If the weight of the best model is low, it is not possible to say that that model really is the best model, meaning other models also have to be considered. In this study, the model with the best AICc and highest weight was used. Significance (p values) for the mixed-effect model was calculated using the package lmerTest (Kuznetsova et al., 2016).

Annual fluxes were calculated using the predict function over the mixed-effect model from library lme4 in R (Bates et al., 2015), and predictions and error terms were back transformed. The predictions were made over half-hourly measurements of soil moisture and soil temperature at 5 cm soil depth in all three microforms just outside the plots.

From these predictions, partitioned fluxes were calculated from the collars without litter as

$$\begin{array}{}\text{(2)}& F_{\mathrm{a}}+\mathit{\epsilon }=(F_{\mathrm{control}}+\mathit{\epsilon })-\left(\right(F_{\mathrm{trench}}+\mathit{\epsilon })-(F_{\mathrm{dead}\mathrm{roots}}+\mathit{\epsilon }\left)\right),\end{array}$$

where F_a is autotrophic CO₂ flux, F_control is the CO₂ flux from the control plots, F_trench is the CO₂ flux from trenched plots, F_dead roots is the CO₂ flux coming from the dead roots in the trenched plots created by the trenching technique and ε is the associated error terms (error terms were propagated).

The annual flux from litter was calculated from the difference in modelled annual fluxes between collars with and without litter.

3 Results

3.1 Temporal trends in soil CO₂ fluxes

Trenching initially led to an increase in soil respiration, followed by a significant reduction in soil CO₂ flux. Soil respiration fluxes from both control and trenched plots showed a clear annual cycle, with the highest fluxes in summer. After the initial perturbation, fluxes from trenched plots are significantly lower than fluxes from control plots (p<0.001), and this difference is greater in the summer (Fig. 3).

Figure 3 Mean CO₂ fluxes from control (grey squares) and trenched (black circles) plots over time, averaging across all microforms (n=12) and both litter treatments. Error bars are standard errors and are in most instances smaller than symbols.

Soil CO₂ fluxes were best explained with a combination of soil moisture, soil temperature, trenching treatment, microform and litter treatment, with an interaction between soil moisture and soil temperature, including plot as a random effect. Table 1 shows model estimates for each variable, with their standard error and p value. The marginal R² was 0.40 and conditional R² was 0.41 on the log scale; since predictions are made in the original units, all uncertainties reinflate. The set of models with a ΔAICc of less than 2 is shown in Table 2. The predictors of the top model were present in over half of the models in the top model set, so this model was used.

Table 1 Model estimates with standard errors and p value for the best-fit model. All numerical predictors were standardized to 1 standard deviation prior to analysis.

Table 2 Model selection summary, showing the four best-ranked models with a delta AICc of less than 2. Models are ranked by AICc and weight, where higher weighted models have more statistical support. Df is the degrees of freedom, LogLik is the log likelihood and ΔAICc is the difference in AICc to the top model.

3.2 Spatial trends and litter contributions to soil CO₂ flux

In both control and trenched plots, fluxes from plough throw (control: 1.23 ± 0.10; trenched: 0.85 ± 0.07 µmol m⁻² s⁻¹; p=0.01) and original surface (control: 1.48 ± 0.10; trenched: 0.83 ± 0.06 µmol m⁻² s⁻¹; p<0.001) were significantly higher than fluxes from furrow (control: 0.90 ± 0.06; trenched: 0.47 ± 0.04 µmol m⁻² s⁻¹). Across all microforms, fluxes from collars with litter (1.03 ± 0.05 µmol m⁻² s⁻¹) were significantly higher than fluxes from collars without litter (0.91 ± 0.05 µmol m⁻² s⁻¹, p=0.008).

3.3 Role of environmental drivers in modulating CO₂ flux

Observed soil CO₂ efflux values correlated positively with soil temperature, whilst soil moisture showed an inconsistent correlation with flux values; a significant (p=0.008) interaction between soil temperature and soil moisture means that at high temperatures CO₂ flux decreases with increasing soil moisture, but at low temperatures flux increases when soil moisture increases, with an inflection between positive and negative moisture dependence between 6 and 7 ^∘C (Fig. 4).

There was no difference in soil temperature between trenched and control plots, but soil moisture was significantly higher in trenched plots than in control plots (p<0.001). To account for the artificially elevated soil moisture in trenched plots, CO₂ fluxes were corrected using the global relationship between soil moisture and CO₂ flux (Fig. 4) to calculate fluxes expected to occur under moisture conditions in control plots.

Figure 4 Combined effect of soil temperature and soil moisture at 5 cm depth on soil CO₂ flux from the control sites, estimated from the top model.

3.4 Partitioned fluxes

Soil CO₂ efflux was partitioned into component fluxes for all measuring dates from August 2014 onwards to remove disturbance-related artefacts observed immediately after trenching. The model prediction of heterotrophic respiration is calculated by subtracting the autotrophic flux from the total soil flux and includes emissions from decomposition of cut roots. Flux simulations based on the soil model details indicate significantly lower autotrophic fluxes than heterotrophic fluxes (p=0.01, Fig. 5). Across all microforms, heterotrophic fluxes represented 61 % and autotrophic fluxes represented 39 % of the total fluxes. From these predictions, annual sums for autotrophic and heterotrophic fluxes have been calculated, giving an average peat decomposition flux of 184 ± 21 g C m⁻² yr⁻¹. Total soil respiration without needle litter is 301.3 ± 34.2 g C m⁻² yr⁻¹, with needle decomposition adding 41.2 ± 53.5 g C m⁻² yr⁻¹ to annual fluxes, giving a total soil respiration including needle litter of 343 ± 35 g C m⁻² yr⁻¹.

Figure 5 Modelled and measured fluxes of heterotrophic (black) and autotrophic (grey) soil CO₂ efflux from the three topographic microforms. Closed symbols are average fluxes with error bars (n=4 plots). Connecting lines are the predicted fluxes using soil temperature and moisture at 5 cm depth. (a) Original surface, (b) plough throw and (c) furrow.

3.5 Impact of litter and roots

Litter fall was 719 ± 71.3 g of litter per square meter per year, with no detectable difference between trenched and control plots. Assuming a C fraction of biomass of 50 % (Mathews, 1993), this represents an input to the soil of 359 g C m⁻² yr⁻¹ via litter fall.

CO₂ flux from surface litter is calculated from the difference in the modelled annual fluxes between the collars with and without litter, with the fluxes from collars with litter significantly higher than the CO₂ flux from the collars without litter (p=0.008, Table 1). C emitted by litter in the control plots appears to be higher than in trenched plots, and the CO₂ flux per gram of litter is 1.7 to 3.6 times higher from the control plots than from the trenched plots (Table 3).

Table 3 Mean amount of C emitted as CO₂ by surface needle litter in litter collars, for both years (standard error in brackets).

Root biomass per square meter with an assumed rooting depth of 25 cm in August 2014 was 761 ± 324, 603 ± 110 and 715 ± 257 g for <2 mm roots, 2–5 mm roots and >5 mm roots, respectively. For both the control and the trench plots, roots smaller than 2 mm declined in total biomass from the start of the experiment to the end of the experiment, but there were no significant differences between the control and trenched plots at the start and end of the experiment. Despite an apparent trend towards higher root biomass in control plots in July 2016, these differences were not statistically significant. There were also no significant differences between the beginning and end of the experiment for root classes 2–5 and >5 mm (Fig. 6).

Figure 6 Root biomass per square meter to 25 cm in control (grey) and trenched (black) plots at the beginning of the experiment (June 2014) and end of experiment (July 2016), split into three root size classes, per microform. (a) Root size <2 mm, (b) root size 2–5 mm and (c) root size > 5 mm. Symbols are average fluxes with error bars representing standard errors (n=4 plots).

3.6 Root decomposition

Root mass in decomposition bags showed a consistent decline over the duration of the experiment (Fig. 7). The calculated decay constant (k) showed systematic differences by microform, with the highest decay rates tending to occur for incubations at the original surface (Table 4).

Figure 7 Root mass remaining in root bags over number of days buried, for each root size class averaged across all microforms. (a) Root size < 2 mm, (b) root size 2–5 mm and (c) root size > 5 mm.

Table 4 Root decay constants and associated C emissions in trenched plots by microform and root size classes (standard error in brackets).

3.7 C flux from dead roots

The amount of C emitted from the decaying roots is calculated using the exponential decay function (1), with the biomass of roots per square meter to a depth of 20 cm in the trenched plots at the beginning of the experiment as M₀. It is assumed that all biomass lost is emitted as CO₂ and that 50 % of roots are C (Mathews, 1993), as conservative assumptions, meaning that estimates are the maximum possible CO₂ flux from decaying roots in trenched plots. To correct for soil CO₂ efflux generated in trenched plots as an artefact of creating dead root biomass, annual estimated CO₂ emissions were corrected by subtracting root-decay-based estimates from trenched plot CO₂ emissions (Table 5).

Table 5 Corrected for dead root decay in trenched plots, heterotrophic soil CO₂ efflux (excluding litter) (F_h) and autotrophic (F_a) fluxes (standard error in brackets) in g C m⁻² yr⁻¹ for both the first year (August 2014–August 2015) and second year (August 2015–August 2016) of the study.

Heterotrophic fluxes represent approximately 24 % and autotrophic fluxes 76 % of the total soil fluxes in the original surface, 37 % and 63 % respectively in the plough throw, 51 % and 49 % respectively in the furrow, and 38 % and 62 % respectively averaged over all microforms.

3.8 Weighted average for soil CO₂ flux in Flow Country forest plantations

In order to scale soil CO₂ fluxes (excluding litter) from different microforms to the level of an entire forest stand, fluxes from individual microforms were weighted by their fractional area (Table 6). This results in a slight shift in proportion of heterotrophic and autotrophic CO₂ flux sources to 40 % and 60 % respectively and a total area-weighted soil CO₂ flux of 289.4 ± 18.6 g C m⁻² yr⁻¹.

Table 6 Area-weighted heterotrophic peat decomposition (F_h) and autotrophic (F_a) fluxes (standard error in brackets) in g C m⁻² yr⁻¹ averaged over both years measured.

4 Discussion

Our results indicate that the soils in these 30-year-old drained and afforested peatlands are a net sink for C, as substantially more C enters the soil as organic matter than is decomposed heterotrophically. The soil surface C balance of the soil under these forest plantations is visualized in Fig. 8, with the annual CO₂ fluxes of the forest plantation based on the area-weighted fluxes. Based on these fluxes, we tentatively estimate the soil C sink (the difference between litter C input and total heterotrophic soil CO₂ efflux) to be 191±64 g C m⁻² yr⁻¹ for our site.

Figure 8 Summary of soil C budget (in g C m⁻² yr⁻¹) under the canopy of a Sitka spruce and lodgepole pine plantation on deep peat. Downward arrows indicate litterfall, whilst upward fluxes show autotrophic (F_a) (black), peat decomposition flux (grey) and surface needle litter (white) across the three main microforms. Total forest soil fluxes with area-weighted values are visualized in the top-right corner, including the estimated soil C sink (grey and white diagonal striped).

We found a C input of 359 g C m⁻² yr⁻¹ via litter fall, similar to other Sitka spruce forests of similar age to these forest plantations in the UK, which range from 273 to 573 g C m⁻² yr⁻¹ (Morison et al., 2012). This is balanced by total soil efflux including litter-derived CO₂ of 342.5 ± 34.7 g C m⁻² yr⁻¹; i.e. the amount of C entering the soil as surface litter alone falls within a similar range to C leaving as CO₂. Missing in the C budget in Fig. 8 is the input of C from root exudates that remain in the soil C pool (net rhizodeposition), potential losses through aquatic C export and root turnover. Gaffney et al. (2020) have measured the aquatic C flux (DOC, DIC, POC) from similar forest plantations in the north of Scotland and have found a flux of 13.9 g C m⁻² yr⁻¹, which is a small flux compared to the CO₂ soil efflux. Root : shoot allocation in forest ecosystems is usually in the order of 1 : 3 (Laiho and Laine, 1997), so it is possible that below-ground litter and exudate deposition could account for a significant C input only partially sampled by our approach (i.e. excluding large root stocks), which means that our tentative soil C sink estimate is conservative.

The rate of peat decomposition in these drained and afforested peatlands is substantial, but overall soil CO₂ efflux and ratio of heterotrophic $/$ autotrophic respiration falls within boreal averages for forested ecosystems (Fig. 9; Bond-Lamberty and Thomson, 2010).

Figure 9 Heterotrophic annual flux against soil respiration annual flux (g C m⁻² yr⁻¹) in boreal forests in black dots and boreal forests on peat soils in grey triangles, from the Bond-Lamberty and Thomson (2010) dynamic database downloaded on 19 April 2020. Our result of this study is added as the red square. See Table A1 in the Appendix for references to studies used to create this graph.

Results from other drained temperate peatland sites under conifer plantation vary from around 260 g C m⁻² yr⁻¹ at a 39-year-old drained Sitka spruce plantation in Ireland on naturally treeless blanket bog (Byrne and Farrell, 2005) to values of 972 g C m⁻² yr⁻¹ on average across 18- to 44-year-old Sitka spruce plantations on poorly drained Dystric Histosols in the south of Ireland (Jovani-Sancho et al., 2018), with similar annual fluxes of 453 g C m⁻² yr⁻¹ reported from under a 43-year-old lodgepole pine plantation in central Scotland (Yamulki et al., 2013). Afforested boreal peatlands also show generally higher soil CO₂ flux values, ranging from 577 to 705 g C m⁻² yr⁻¹ (e.g. Mäkiranta et al., 2008; Minkkinen et al., 2018; Ojanen et al., 2010). Also taking into account possible measurement bias associated with the chemical absorption method used by Byrne and Farrell (2005), which generally underestimates fluxes compared to more reliable infrared gas analysis based methods (Janssens and Ceulemans, 1998), our results of 342.5 g C m⁻² yr⁻¹ sit at the low end of total soil respiration fluxes. These low rates are corroborated by independent chamber flux measurements across forest sites in the Flow Country (Hermans et al., unpublished data) and are likely to relate to poor site fertility and specific biochemical properties of the peat at these sites (Hermans et al., 2018).

Our peat oxidation rates of 115 ± 16 g C m⁻² yr⁻¹ are also low compared to temperate afforested peatlands, which have been reported to vary between 211 and 407 g C m⁻² yr⁻¹ for Sitka spruce sites in Ireland (Jovani-Sancho et al., 2018), whilst boreal peatlands show an even larger range of heterotrophic flux sums (145 to 670 g C m⁻² yr⁻¹; Ojanen et al., 2010). Our low estimates result in part from the root decomposition correction, which reduced heterotrophic flux estimates by around 60 g C m⁻² yr⁻¹ and is not routinely applied across other studies.

Minkkinen et al. (2018) measured carbon fluxes of a drained, naturally forested peatland (Kalevansuo peatland forest) in the south of Finland over 4 years. They found that the ecosystem was a strong sink of CO₂ in all studied years, with an average net ecosystem exchange (NEE) for the 4 years of −234 g C m⁻² yr⁻¹. By subtracting the carbon sink of the tree stand from eddy-covariance-derived NEE the authors show that the soil was a carbon sink of −60 g C m⁻² yr⁻¹. However, modelling their forest soil respiration fluxes from chamber measurements, they found the peat respiration alone made up 53 % of the total forest floor respiration flux, litter 22 %, roots 16 % and autotrophic respiration of above-ground vegetation 8 %, with soils representing a weak C source. Our results show a lower percentage of peat respiration of 38 % of the total soil respiration minus the litter flux. The higher percentage found by Minkkinen et al. (2018) could possibly be explained by them not applying a dead root correction to their fluxes, which as we show could lead to a big difference in fluxes. For forestry-drained boreal peatlands, Ojanen et al. (2013) report a range of soil C balances from weak source to sinks depending on site fertility. Accounting for below-ground litter production and fine root turnover is an important source of uncertainty in soil C balance estimates (Ojanen et al., 2014), and whilst our soil C balance focuses on surface litter only, our calculations indicate a weak C sink even before below-ground litter input is considered.

Trenching is likely to underestimate heterotrophic respiration and rates of peat mineralization since this approach does not account for rhizosphere effects such as priming (Walker et al., 2016). Therefore, it is likely that the estimated rate of peat oxidation from trenching is conservative. The observed difference of 1.7 to 3.6 times higher CO₂ flux per gram of litter from the control plots than from the trenched plots (Table 3) indicates that heterotrophic processes are reduced under trenching. In presence of an active rhizosphere (control plots), decomposition of needle litter and/or soil organic matter (SOM) appears to be faster than when the rhizosphere is not active (trenched plots). This priming of organic matter decomposition is likely to be the result of microbial activity in the soil stimulated by the addition of easily accessible C from roots (Kuzyakov, 2006). Friggens et al. (2020) recently demonstrated the accelerated loss of organic soil under coniferous and deciduous forest plantations on peaty podzols in Scotland, which is consistent with soil priming through the elevated input of C from trees in the rhizosphere. Several studies have shown that mycorrhizal fungi are involved in soil priming (Kuyper, 2017; Paterson et al., 2016), further underlining the interaction between plant-associated biotic activity and decomposition in underlying peat. Therefore, in the control plots a slightly larger proportion of the total CO₂ flux could be heterotrophic decomposition than the fluxes from the trenched plots suggest, which means there could be an underestimation of heterotrophic flux in our results, in line with results from the literature (Subke et al., 2004, 2011).

To calculate the root biomass at the start of the experiment, one soil core per microform was taken, and since trees were standing close to each other (1.5 to 2 m apart) this was assumed to be representative for the whole microform. It was not possible to distinguish between living and dead roots in the soil cores, so our results might overestimate the living root biomass. The dead root emission correction made a big difference to the ratio of heterotrophic to autotrophic flux, going from 61 % and 39 % respectively over all microforms to 38 % and 62 % respectively (without area-weighting of fluxes), so a decrease in heterotrophic flux of 23 %. This is in line with the corrections used in other studies; a comparison of corrections used in trenching studies indicates a range from 2 % to 24 %, with an average of 12 % (Subke et al., 2006). This big difference in the fraction heterotrophic : autotrophic flux suggests that even 2 years after trenching, decaying root biomass makes significant contributions to the CO₂ flux.

That soils under a 30-year-old conifer plantation on deep peat are a C sink is an interesting finding. However, in order to determine if these soils are a long-term sink or source of carbon, the balance between soil C input from roots and litter and C loss via oxidation should be quantified over the lifespan of the plantation. The average peat depth in these forest plots is 126 (±16) cm, with 0.47 kg C m⁻² per centimetre depth (Cannell et al., 1993), which represents ca. 59.3 (±7.3) kg C m⁻² stored in the peat under these plantations. The largest carbon losses most likely occurred in the initial planting phase of afforestation (Simola et al., 2012; Vanguelova et al., 2019; Zerva and Mencuccini, 2005), but fluxes have not been measured through this phase and cannot be quantified as new policy prevents planting on deep peat (Forestry Commission, 2016). In other parts of the world, peatland drainage is still actively happening, and studies from these sites show very high rates of peat oxidation during the first 5–10 years of conversion (e.g. Cooper et al., 2020; McCalmont et al., 2021; Minkkinen et al., 2018; Prananto et al., 2020). In the UK, forest plantations on deep peat sometimes end in clear felling of the site and restoration of the peat. The results of this study could also help inform what the CO₂ fluxes will be when timber is harvested and large quantities of felling residues are left above ground as well as whole root systems below ground. However, we note that changes in water table and soil moisture conditions created by the soil rewetting associated with clear felling will have significant and separate impacts beyond the conditions of active drainage under which we took our measurements.

Appendix A

Table A1 List of articles used to create Fig. 9.