Carbon Budget of Undrained and Drained Nutrient-Rich Organic Forest Soil

Abstract

The impact of the moisture regime on the carbon budget of organic soils with different nutrient statuses has not been fully studied in hemiboreal forests thus far. This study evaluated soil carbon (C) stock changes in forests with drained and undrained nutrient-rich organic soils by estimating C loss through respiration and C input through the litter. The study sites included forest stands dominated by Norway spruce (Picea abies), silver birch (Betula pendula), black alder (Alnus glutinosa), and clear-cuts. Soil respiration was measured using the chamber method, and to estimate the soil C input by litter—the biomass and the C content of the foliar litter, ground vegetation, and fine-root production were measured. The soil in forest stands acted as a C sink. The carbon dioxide (CO₂) removal rates of 0.4 ± 0.4 t C ha⁻¹ year⁻¹ and 0.1 ± 0.4 t C ha⁻¹ year⁻¹ were estimated for undrained and drained soil in forest stands, respectively. The soil in the clear-cuts acted as a CO₂ source, and the annual emissions ranged from 0.4 ± 0.4 t C ha⁻¹ year⁻¹ in undrained to 0.9 ± 0.7 t C ha⁻¹ year⁻¹ in drained conditions. The reason for the soil in clear-cuts being a C source was increased C loss by respiration and reduced soil C input by litter. Furthermore, the mean soil C input by ground vegetation biomass in the clear-cuts was considerably higher than in the forest stands, which did not compensate for the increase in soil respiration and the absence of C input by foliar litter and the fine roots of trees. The results of the study on annual soil C stock changes can be used as an emission factor in national greenhouse gas inventories of forest land in the hemiboreal zone.

1. Introduction

In accordance with the Paris Agreement, the European Union, including Latvia, has committed to achieving climate neutrality by 2050. These policy targets promote seeking forest management practices that contribute to C sequestration [1,2,3,4], reflected in the topicality of related studies, including the estimation of ecosystem greenhouse gas (GHG) balance. Reducing the GHG emissions from the main national sources, such as the transport, energy, and agriculture sectors, which currently account for around 88% of Latvia’s total GHG emissions, will not be sufficient to achieve the climate neutrality target set by Paris Agreement. To compensate for the irreducible GHG emissions in these sectors, the land use, land-use change, and forestry (LULUCF) sector must ensure the equivalent sequestration of carbon dioxide (CO₂). Boreal forests are often identified as ecosystems with carbon (C) sequestration potential [5]. Therefore, forest land is the only land-use category of the LULUCF sector in which an increased rate of CO₂ sequestration by the implementation of climate change mitigation measures has the potential to offset the country’s total GHG emissions.

Efforts to achieve the C sequestration and GHG mitigation potential of the forest ecosystem can be implemented with targeted activities that promote C sequestration in biomass, soil, and harvested wood products, as well as by replacing fossil fuels with biomass. The role of forest management in tree biomass C sequestration is well understood and modeled, but understanding the process of C sequestration and assessing changes in the C stock of the soil remain limited [6]. Previous studies show that climate change mitigation measures targeted at organic soil management are often the most effective, but the climate change mitigation potential of organic soils is not being used⁠ fully [7,8]; this is largely related to the lack of knowledge. The most frequently identified climate change mitigation measures related to the management of organic forest soils are the afforestation of organic soils and the restoration of the natural moisture regime [9]. However, there is a lack of scientific evidence that the restoration of the natural moisture regime of organic soils promotes reductions in GHG emissions and an increase in the C sink of the boreal forest ecosystem in the current climate conditions. GHG emissions from undrained forest soil are not within the scope of the national GHG inventory (Inventory) reports, hindering interest in such studies. However, quantitative awareness of such emissions is crucial for the comparison of carbon stock change (CSC) of drained organic soil relative to undrained soil to fully understand the climate impact of drainage and to enable possibilities of implementing the potentially most effective climate change mitigation measures in forest land management. The drainage of organic soils is often considered a climate-harmful management practice, although knowledge of annual soil GHG emissions is highly uncertain [8]. Currently, there is a lack of common understanding of the impact of soil drainage on forest ecosystem GHG emissions and the C balance. Some studies indicate that the drainage of nutrient-poor organic soils in boreal forests has a significant impact on ecosystem CO₂ sequestration [10], while the drainage of nutrient-rich organic soils may turn forest ecosystems into GHG emission sources when soil C and nitrogen (N) loss is not compensated by increased forest growth [11].

The most commonly mentioned shortcoming of previous scientific articles on net CO₂ emissions from forests with drained organic soils in boreal and temperate climate regions is a necessity to subtract below- and above-ground biomass respiration from the reported results and incorporate litter production or decomposition rates [8]. Thus, the results reported require further processing or additional data to enable the quantification of annual soil CSC. Another shortcoming is the uneven site spatial coverage of the previous studies. Most of the organic soil CSC estimate results were obtained from drained boreal peatland studies carried out in Finland, while most of the study sites representing a temperate zone are located in the southern part of Sweden [8]. The results of organic soil CSC estimates in the Baltic states representing hemiboreal forests are reported by four articles on drained peatlands [12,13,14,15]. Despite the fact that the availability of study results on drained organic forest soils has increased, they are still scarce, considering the variability of the factors affecting CSC in forest ecosystems. No CSC estimates of undrained organic forest soils have been reported in the region, as studies on undrained organic soils are usually carried out in pristine or recently recultivated peatlands.

According to the acknowledgment that there is a lack of studies evaluating the impact of different long-term soil moisture regimes on soil CSC [8] and observations that organic forest soil CO₂ emissions can be comparably higher in forest sites with increased soil fertility [16], room for improvement in the Inventory and capabilities to plan climate change mitigation measures is recognized. The estimated GHG emissions of 1.7 million t CO₂ equivalents (14.4% of the total emissions of Latvia) from drained organic soil in the forest land category in 2020 [17] show the significance of accurate organic soil emission estimates in the national Inventory. The currently applied country-specific emission factor (0.52 t C ha⁻¹ year⁻¹) for the estimation of CO₂ emissions from drained organic soil in forest lands in Latvia is developed by the C stock inventory method conducted in forests with nutrient-poor to moderately rich (Callunosa turf. Mel., Vacciniosa turf. Mel. Additionally, Myrtillosa turf. Mel.) soils according to the national forest site type classification [18]. The country-specific emission factor is applied to all organic soils in forests, while according to the national forest inventory, the share of drained (17%) and undrained (4%) forest site types with nutrient-rich organic soils, where potentially higher soil CO₂ emissions may be expected, is 21%. Therefore, the currently used emission factor may introduce accuracy errors in the estimations.

This study aimed to estimate the CSC of drained and undrained nutrient-rich organic soils using empirical data on soil CO₂ emissions and soil C input by:

• Foliar litter (LF);
• Ground vegetation (above- and below-ground biomass of herbs and grasses, GV);
• Fine roots of trees (FR);
• Moss and dwarf shrubs.

This study contributes to the improvement of the national GHG inventory and provides a scientifically valid assessment of potential soil drainage effects on CO₂ emissions to support decision-making on climate change mitigation measures.

2. Materials and Methods

2.1. Study Site Description

The study was carried out in central Latvia (Figure 1) on the forest stands of a hemiboreal zone with undrained (Dryopterioso-caricosa and Filipendulosa) and drained (Oxalidosa turf. mel.) forest site types characterized by nutrient-rich organic soil. For the forest stands to be accepted as study sites, the compliance with drainage status, the average peat layer depth (>30 cm in undrained sites and >20 cm in drained sites), and the characteristic vegetation, as defined in the national forest site type classification (

Table 1. Dominant ground vegetation in the study sites.

Forest Site Type	Ground Vegetation
Dryopterioso-caricosa	Thelypteris palustris Schott, Carex (L.), Iris pseudacorus (L.) Fuss, Scirpus (L.), Lysimachia vulgaris (L.), Cirsium oleraceum (L.) Scopoli, Filipendula ulmaria (L.) Maximowicz, Angelica sylvestris (L.)
Filipendulosa	Filipendula ulmaria (L.), Urtica dioica (L.), Geum rivale (L.), Paris quadrifolia (L.), Caltha palustris (L.), Solanum dulcamara (L.)
Oxalidosa turf. mel.	Cirsium oleraceum (L.) Scopoli, Hylocomium splendens (H.) Schimper, Rhytidiadelphus, Brachythecium, Vaccinium myrtillus (L.), Dryopteris filix-mas (L.) Schott

), were analyzed. One round sample plot (500 m²) was established in each of the selected study sites. The distance to the nearest drainage ditches from the sample plots was at least 300 m and 100 m in the study sites with undrained and drained soil, respectively.

During the collection of the empirical data (from October 2019 to June 2021), the air temperature in the study sites ranged from 8.0 ± 0.7 °C to 31.4 ± 0.1 °C (mean 9.2 ± 0.8 °C) and the annual precipitation ranged from 472 mm to 860 mm (average 668 ± 136 mm) according to the data provided by the meteorological stations of the Latvian Environment, Geology, and Meteorology Centre (distance from study site less than 30 km).

The scope of the study included 31 forest stands in total, with the dominant tree species being Norway spruce (Picea abies (L.) Karsten), silver birch (Betula pendula Roth), and black alder (Alnus glutinosai (L.) Gärtner) at different stages of stand development (hereafter spruce, birch, and alder, respectively), from clear-cuts to mature stands (

Table 2. Characteristics of forests sites.

Parameter	Undrained Forest Sites				Drained Forest Sites
	Spruce	Birch	Alder	Clearcut	Spruce	Birch	Alder	Clearcut
Number of study sites	1	3	5	1	12	3	2	4
Age, years	67	21–77	10–80		14–86	18–60	26–53
Diameter, cm	31	12–29	4–23		2–27	9–27	17–24
Height, m	28	12–28	4–29		2–24	9–22	17–26
Basal area, m2 ha−1	61	17–71	8–57		8–72	19–60	32–56
Growing stock, m3 ha−1	335	78–365	35–325		7–521	38–210	123–254
Thickness of peat layer, cm	68	31–52	30–99	47	37–99	25–75	60–70	63–99

The table shows the range of characteristics of forest sites.

). The results of the individual peat-layer thickness measurement replicates varied from 23 cm to ≥100 cm (mean 75 ± 7 cm) in the undrained sites and from 25 cm to ≥100 cm (mean 54 ± 12 cm) in the drained sites.

2.2. Sampling and Laboratory Analysis of Soil and Soil Flux

Soil CO₂ flux monitoring was conducted using the manual closed static nontransparent chamber method [19] for 12 consecutive months. Chamber collars were installed at a depth of 5 cm in 5 replicates in each study site. During collar installation, root damage and disturbance of the litter layer were avoided as much as feasible, and GV was left intact throughout the whole monitoring period. Therefore, the flux monitoring represents the CO₂ exchange between the soil surface (including vegetation enclosed in the chamber) and the atmosphere, the sum of soil heterotrophic respiration, and the autotrophic respiration of roots and aboveground ground vegetation (Rfloor), respectively. Soil flux was sampled with an interval of 4 weeks from the chambers in each of the collar positions immediately and at 10, 20, and 30 min after positioning the chambers on the collars. The samples were collected using underpressurized (0.3 mbar) glass vials and tested using the gas chromatography method [20]. The atmosphere and soil temperature at a 5 cm depth (Ts), as well as the groundwater level (using a PVC pipe installed up to a depth of 140 cm), were recorded during the soil flux sampling.

The soil samples were collected with 100 cm³ cores from fixed soil depths of 0–10 and 10–20 cm in 2 replicates [21]. The soil samples were prepared according to LVS ISO 11464:2005, and the bulk density was determined according to LVS ISO 11272:2017. The soil chemical parameters were determined using standard methods (

Table 3. Standard methods used in analyzing soil samples.

Parameter	Unit	Method Principle	Standard Method
Bulk density	kg m−3	Gravimetry	LVS ISO 11272:2017
Total C	g kg−1	Elementary analysis (dry combustion)	LVS ISO 10694:2006
Total N	g kg−1	Elementary analysis (dry combustion)	LVS ISO 13878:1998
Carbonate (CaCO3)	g kg−1	Volumetry	LVS EN ISO 10693:2014
pH	unit	Potentiometry	LVS ISO 10390:2021
HNO3 extractable potassium (K), calcium (Ca), magnesium (Mg) and phosphorus (P)	g kg−1	ICP-OES	LVS EN ISO 11885:2009

). The content of organic C was calculated by subtracting the value of carbonate C from the total C value. In addition, the soil organic C/total N ratio (C/N ratio) was calculated as a proxy to characterize the decomposition of soil organic matter.

2.3. Estimation of Soil Respiration

The acquired analysis results of the CO₂ concentration in the chambers during soil flux sampling were used to calculate the slope values of the linear regression equations characterizing the gas concentration changes over time. The instantaneous Rfloor was calculated using the following equation:

$$Rfloor=\frac{M P V slope}{R T A}$$

(1)

where Rfloor is the instantaneous Rfloor, µg CO₂ m² h⁻¹; M is the molar mass of CO₂, 44.01 g mol⁻¹; R is the universal gas constant, 8.314 m³ Pa K⁻¹⋅mol⁻¹; P is the assumption of air pressure inside the chamber, 101,300, Pa; T is the air temperature, K; V is the chamber volume, 0.063 m³; slope is the CO₂ concentration changes over time, ppm h⁻¹; and A is the collar area, 0.1995 m².

The annual Rfloor was estimated by summing the calculated hourly Rfloor. We calculated the hourly Rfloor by interpolating the measured instantaneous Rfloor using the R₁₀ and Q₁₀ parameters [22,23,24], the relationship between the atmospheric temperature and Ts evaluated within the study and the hourly average air temperature data from the nearest meteorological stations. The hourly Rfloor was calculated using the following equation:

$$Rfloor=R_{10}{Q}_{10}^{\frac{T_S-10}{10}}$$

(2)

where Rfloor is the hourly Rfloor, kg CO₂ ha⁻¹ h⁻¹; R₁₀ is the Rfloor at a soil temperature of 10 °C, kg CO ha⁻¹ h⁻¹; Q₁₀ is the temperature sensitivity; and T_S is the soil temperature, °C.

The following equation (R² = 0.81, p < 0.001), elaborated by the results from previous studies [25], was used to recalculate the annual Rfloor to soil heterotrophic respiration (Rhet):

$$\ln \left(Rhet\right)=1.22+0.73\ln \left(\mathrm{Rs}\right)$$

(3)

where Rhet is soil heterotrophic respiration, t CO₂ ha⁻¹ year⁻¹, and Rs is soil respiration, t CO₂ ha⁻¹ year⁻¹.

The annual Rfloor and Rhet were estimated by stratifying the empirical data acquired in the study according to soil moisture regime (undrained and drained), forest land status (forest stand or clear-cut), and forest type (deciduous or coniferous) to allow the application of study the results for the improvement of the national GHG inventory.

2.4. Estimation of Soil C Input by Litter

The LF samples for the estimation of the annual LF biomass were collected using five conically shaped litter traps (surface area 0.5 m²) installed in each study site according to the manual methods and criteria for harmonized sampling, assessment, monitoring, and analysis of the effects of air pollution on forests, prepared on behalf of the Programme Co-ordinating Centre and Task Force of ICP Forests [26]. The samples were collected for 12 consecutive months with an interval of 4 weeks.

Separate above- (aGV) and below-ground ground vegetation (bGV) samples were collected in 4 replicates from 20 cm × 20 cm square fields in each study site at the end of the vegetation season when the vegetation biomass had peaked [23]. The FR production samples were collected using the modified ingrowth core method based on a flexible polyester cylindrical bag (diameter 35 mm) with a mesh size of 2 cm × 2 cm installed 60 cm deep in the soil in three replicates in each study site [27,28]. Mesh bags were installed in autumn and removed from the soil after a year by cutting the roots around the bag. The roots of trees were removed from the collected bGV samples, while the roots of GV were removed from the collected FR samples. The soil particles from both types of below-ground biomass samples were removed by wet sieving.

The litter sample dry matter was determined by oven drying (70 °C) the samples; the C content was analyzed by dry combustion using an element analyzer according to LVS ISO 10694:2006. It was assumed that the biomass of the collected foliar litter, GV, and FR production was equal to the annual mortality and respective soil C input:

$$\mathrm{Cinput}=\frac{\mathrm{m}\times 10000\times \mathrm{C}}{\mathrm{S}\times 100}$$

(4)

where Cinput is the annual soil C input by litter, t ha⁻¹ year⁻¹; m is the dry matter of litter, t; C is the C content of litter, %; and S is the area of the litter sampler (cross-sectional area of root ingrowth bag, area of LF trap, and area of GV collection field), m².

2.5. Estimation of Forest Soil Annual CSC

The soil CSC was calculated as the sum of the soil C input by annual biomass mortality (LF, GV, FR, mosses, and shrubs) and soil C loss by Rhet. The estimated soil CSC was expressed as the mean annual CSC within 240 years of forest management in a business-as-usual scenario. Assumptions of yearly stand age and basal area development within a period of 240 years of forest management (including the impact of harvesting), which were used as variables for the annual soil CSC calculations (Figure 2), are based on the National Forest Inventory data and national stand growth models [29,30,31,32].

The assumptions of the yearly dynamics of the basal area were used as variables for the calculation of the annual soil C input according to the study results. The annual LF C input was calculated using equations for the relationship between the basal area and C content in the annual LF biomass, while the annual soil C input by FR and GV was estimated according to forest land status (forest stand or clear-cut). It was assumed that forest stand or clear-cut status could be determined by the national stand basal area thresholds identifying unproductive stands: 6, 4, and 5 m² ha⁻¹ for spruce, birch, and alder stands, respectively.

The study results for soil C input and Rhet are supplemented by data on the annual litter biomass of dwarf shrubs and mosses. The biomass of shrubs and mosses was calculated using the equations elaborated within a study conducted in boreal forests [33]:

$$B{\left(spruce\right)}_{shrubs}=\left(10.375-0.033A+0.001A^2-0.000004A^3{)}^2-0.5\right)$$

(5)

$$B{\left(broadleves\right)}_{shrubs}=\left(7.102-0.0004A^2{)}^2-0.5\right)$$

(6)

$$B{\left(spruce\right)}_{mos}=\left(19,282+0.164A-0.000001A^3{)}^2-0.5\right)$$

(7)

$$B{\left(broadeleves\right)}_{mos}=\left(13.555+0.056A{)}^2-0.5\right)$$

(8)

where B(spruce)_shrubs and B(broadleaves)_shrubs are the aboveground biomass of shrubs in coniferous stands and broadleaves forests (kg ha⁻¹), respectively. B(spruce)_moss is the aboveground biomass of moss in coniferous forests (kg ha⁻¹), and B(broadleaves)_moss-is the aboveground moss biomass in deciduous forests (kg ha⁻¹). A is stand age (years).

The annual soil C input by shrubs and mosses was calculated with the assumption that the share of C in the biomass was 47.5% [34] and by multiplying the biomass values with a turnover rate of 0.25 and 0.33 for shrubs and mosses, respectively [33]. It was assumed that 70% of the total C input by dwarf shrubs and mosses contributed to belowground biomass mortality [35,36,37].

2.6. Statistical Analysis

Statistical analyses were carried out using R (R version 4.0.3; RStudio version 2022.07.1 + 554). A Mann–Whitney U test was used to compare the differences between the two data groups. The correlations were tested with Spearman’s ρ. A significance level of α = 0.05 was applied in all the tests. The uncertainty of the study results was expressed with confidence intervals (α = 0.05).

3. Results

3.1. Soil Characteristics of the Study Sites

The mean organic C content in the top 20 cm of the soil in the studied stands with drained soil was 48.7 ± 4.0% and 45.5 ± 4.3% in stands with undrained soil. Thus, the soil in the studied stands complies with the definition of organic soil [38]. The mean soil C/N ratios in the drained and undrained soil were 19.4 ± 2.8 and 19.2 ± 2.9, and the soil bulk densities were 420 ± 40 kg m⁻³ and 435 ± 43 kg m⁻³, respectively.

3.2. Soil Respiration

During the soil CO₂ flux monitoring period, the measured Ts ranged from −1.3 to 22.3 °C, while the instantaneous Rfloor ranged from 0.6 to 97.8 µg C m⁻² s⁻¹ (Figure 3). The highest mean instantaneous emissions were found in clear-cuts. The difference between the measured mean Rfloor in the clear-cuts with drained (31.5 ± 7.0 µg C m⁻² s⁻¹) and undrained (33.4 ± 14.4 µg C m⁻² s⁻¹) soil is not significantly different. The measured mean Rfloor in the forest stands with different dominant tree species, and the soil moisture regimes were also not significantly different from each other and ranged from 18.5 ± 7.1 µg C m⁻² s⁻¹ in spruce stands with undrained soil to 25.9 ± 7.2 µg C m⁻² s⁻¹ in birch stands with drained soil. However, the difference between the measured mean Rfloor in the clear-cuts (31.9 ± 62 µg C m⁻² s⁻¹) and forest stands (21.7 ± 1.8 µ C m⁻² s⁻¹) is significantly different.

The relationship between the Ts and Rfloor can be expressed by exponential regression. The equation (Rfloor = a $\times$ exp (b $\times$ Ts) coefficients a and b are summarized in

Table 4. Summary of the models for the prediction of soil Rfloor by soil temperature at a depth of 5 cm.

Moisture Regime	Forest Site Characteristics	Characteristics of Measured Data (Range)		Model Coefficients		Model Characteristics
		Ts, °C	Rfloor, µg C m−2 s−1	a	b	RMSE	Q10
Drained	Clear-cut	−1.3…19.0	0.5…77.2	5.784	0.141	11.9	4.10
	Deciduous	0.5…22.3	1.2…97.8	4.476	0.143	10.3	4.18
	Coniferous	−0.6…18.2	4.2…59.7	6.235	0.118	6.2	3.25
Undrained	Clear-cut	1.1…21.3	6.3…78.5	7.298	0.124	10.6	3.46
	Deciduous	0.7…19.3	0.6…61.3	4.700	0.136	8.5	3.90
	Coniferous	1.4…17.4	6.2…45.9	5.798	0.124	3.0	3.46

. According to the study results, the Rfloor tended to be more sensitive to Ts changes in deciduous forests and drained clear-cuts (Q₁₀ from 3.9 to 4.18) as compared to the other study sites, namely, clear-cuts with undrained soils and coniferous forests (Q₁₀ from 3.25 to 3.46).

The annual Rfloor was estimated by applying the prediction models developed by the study and hourly air temperature data within the study period in combination with the observed relationship between the temperature of the atmosphere and Ts (R² = 0.87, p < 0.05):

$$\mathrm{Ts}=0.715t_{air}+1.719$$

(9)

where Ts is the soil temperature at a 5 cm depth, °C, and t_air is the air temperature, °C. The estimated annual Rfloor ranged from an average of 5.1 ± 0.2 t C ha⁻¹ year⁻¹ and 5.1 ± 2.6 t C ha⁻¹ year⁻¹ in alder stands with drained soil and birch stands with undrained soils, respectively, to 7.9 ± 3.3 t C ha⁻¹ year⁻¹ in clear-cut stands with undrained soil. The estimated annual mean Rfloor in the forest sites and clear-cuts was 6.2 ± 0.4 t C ha⁻¹ year⁻¹ and 7.7 ± 1.7 t C ha⁻¹ year⁻¹, respectively. The empirical data acquired show a correlation between the Rfloor and GV biomass (r was 0.4 to 0.55 for bGV and GV, respectively). It was also observed that the soil C content (r = 0.51) and LF biomass (r = −0.59) had a moderate correlation with Rfloor, while stand age had a weak (r = −0.36) but significant (p < 0.05) impact on Rfloor (Figure 4).

3.3. Soil C Input by Litter

In the study sites, the aGV (mean C content 47.1 ± 0.7%) and bGV (mean C content 49.7 ± 0.8%) at the end of vegetation season ranged from 6.1 to 8.2 t ha⁻¹ (average 6.9 ± 1.0 t ha⁻¹) in clear-cuts to 1.3 to 6.5 t ha⁻¹ (average 3.5 ± 0.7 t ha⁻¹) in forest stands. While there was no statistically significant relationship identified between the aGV and soil chemical parameter data, the bGV data had a relationship with the parameters indicating soil fertility. The bGV data had a moderate correlation with soil N (r = 0.51), K (r = 0.49), P (r = 0.69) content, and C/N ratio (r = −0.62). Although GV had a moderate correlation (r = 0.51) with stand age, due to the lack of study data available to elaborate models based on stand variables, soil CSC modeling was chosen to be performed by fixed ground vegetation biomass values stratified according to forest land status (forests stand or clear-cut), moisture regime (drained or undrained soil), and dominant tree species (coniferous or deciduous).

The same approach was applied regarding the FR litter data. The estimated FR production in the forest stands ranged from 0.1 to 1.8 t ha⁻¹ (average 0.8 ± 0.2 t ha⁻¹). Although moderate correlations between the estimated FR production and soil fertility characteristics data exist, these relationships were not found to be significant, except in the case of soil Mg content. The study data show that a lower annual average groundwater level tended to increase FR production (r = 0.38).

The estimated annual LF biomass ranged from 0.5 to 5.7 t ha⁻¹ (average 3.3 ± 0.5 t ha⁻¹) with a mean C content of 52.1 ± 0.2%. The annual foliar litter biomass data had a moderate (r from 0.44 to 0.65) correlation with average tree diameter, basal area, height, and growing stock (in order of increasing correlation) to a high correlation with stand age (r = 0.84). The basal area was chosen as a predictor for the explanation of the annual foliar litter biomass due to its better representation of the impact of deciduous or coniferous tree species. The annual soil C input by litter in the spruce stands had a linear relationship with the stand basal area; in the study sites, the estimated C input increased from 0.26 to 2.34 t C ha⁻¹ year⁻¹ in stands with a basal area of 8 m² ha⁻¹ and 45 m² ha⁻¹. The acquired data suggest that, in the case of deciduous forests, the LF stands had a steeper biomass increase until the basal area reached around 20 m² ha⁻¹. When the annual C input by litter reaches around 1.5 t C ha⁻¹ year⁻¹, further increases in the basal area have a more gradual impact on litter biomass increases (Figure 5).

3.4. Summary of Estimated Annual Soil CSC

According to the empirical data of the annual average Rfloor recalculated to Rhet and the soil C input estimated in the study, summarized in

Table 5. Summary of estimated soil respiration and soil C input (t C ha⁻¹ year⁻¹) in the study sites.

Moisture Regime	Forest Site Characteristics	Rfloor	Rhet	aGV	bGV	FR	LF	Net Balance
Drained	Clear-cut	−7.6 ± 2.2	−4.3 ± 0.90	2.2 ± 0.3	1 ± 0.2		0.2 ± 0.2	−0.9 ± 0.7
	Deciduous	−6.2 ± 1.2	−3.7 ± 0.53	0.7 ± 0.3	1.3 ± 0.3	0.3 ± 0.3	1.6 ± 0.3	0.2 ± 0.4
	Coniferous	−6.3 ± 0.5	−3.7 ± 0.23	0.8 ± 0.2	0.8 ± 0.2	0.4 ± 0.1	1.7 ± 0.5	0.0 ± 0.3
Undrained	Clear-cut	−7.9 ± 1.2	−4.4 ± 0.50	2.4 ± 0.3	1.4 ± 0.3		0.1 ± 0.4	−0.4 ± 0.4
	Deciduous	−6.1 ± 0.7	−3.7 ± 0.32	0.7 ± 0.3	1.0 ± 0.4	0.3 ± 0.1	1.7 ± 0.6	0.0 ± 0.4
	Coniferous	−5.1 ± 1.2	−3.3 ± 0.50	0.5 ± 0.3	0.8 ± 0.3	0.5 ± 0.1	2.2 ± 0.4	0.8 ± 0.4

, the soil C loss by Rhet is compensated by the annual soil C input in forest stands with drained and undrained soils, while drained and undrained soil in clear-cuts is a net CO₂ source.

According to the modeling exercise explained in Section 2.4., within a 240-year forest management cycle, the annual soil CSC ranged from −1.0 to 2.6 t C ha⁻¹ year⁻¹ (mean 0.4 t C ha⁻¹ year⁻¹) in deciduous forests and from −0.6 to 2.9 t C ha⁻¹ year⁻¹ (mean 0.7 t C ha⁻¹ year⁻¹) in coniferous forests with undrained soil, whereas in forests with drained soil, the annual net C balance ranged from −1.3 to 1.0 t C ha⁻¹ year⁻¹ (mean 0.1 t C ha⁻¹ year⁻¹) in stands with deciduous-dominant species and from −1.3 to 0.1 t C ha⁻¹ year⁻¹ (mean −0.6 t C ha⁻¹ year⁻¹) in stands with coniferous-dominant species (Figure 6). These results indicate that long-term drainage reduces C uptake by nutrient-rich organic soil in managed forests. After drainage, the soil in deciduous forests may remain C-neutral, but in coniferous forests, the soil may become a CO₂ source.

The soil CSC modeling results showed that aboveground litter, on average, contributed to the annual soil C input by 60 ± 4%, of which 40 ± 10% was from LF and 57 ± 9% from aboveground ground vegetation. The main component of the soil C input was GV, contributing 60 ± 7% of the total annual soil C input (

Table 6. Summary of estimated forest management cycle’s annual average soil CSC (t C ha⁻¹ year⁻¹) for the study site measurements.

Soil CSC Component	Drained		Undrained
	Deciduous	Coniferous	Deciduous	Coniferous
Rhet	−3.82 ± 0.45	−3.85 ± 0.36	−3.88 ± 0.7	−3.54 ± 1.38
LF	1.14 ± 0.21	0.81 ± 0.24	1.12 ± 0.46	0.73 ± 0.22
FR	0.31 ± 0.25	0.44 ± 0.13	0.32 ± 0.29	0.53 ± 0.16
aGV	0.93 ± 0.61	1.1 ± 0.31	1.79 ± 0.61	1.59 ± 0.45
Aboveground shrubs	0.01 ± 0.01	0.01 ± 0.01	0.01 ± 0.01	0.01 ± 0.01
Aboveground mosses	0.02 ± 0.02	0.09 ± 0.04	0.02 ± 0.02	0.11 ± 0.04
bGV	1.27 ± 0.42	0.59 ± 0.17	1.13 ± 0.52	1.00 ± 0.28
Belowground shrubs, mosses	0.07 ± 0.07	0.25 ± 0.1	0.07 ± 0.07	0.28 ± 0.12
Net balance	−0.08 ± 0.03	−0.55 ± −0.12	0.43 ± 0.17	0.71 ± 0.25

).

4. Discussion

4.1. Method of Rhet Calculation

Previous studies in boreal forests mainly focused on the direct evaluation of Rhet [8]. In our study, the decision to estimate soil respiration by Rfloor measurements is a result of a methodological compromise allowing the acquisition of collected gas sample analysis results for soil CH₄ and N₂O flux estimates as well [39]. For this reason, in our study, the Rfloor measurements were recalculated to Rhet using an equation elaborated from the data of previous studies. While such an approach may introduce additional uncertainty in the results of Rhet estimates, it allowed us to maintain the simplicity of the soil CSC calculation. Although the use of direct Rhet measurement results would avoid such unknown potential errors, the method (trenching) required to acquire such Rhet data may introduce other errors in the flux estimates due to altered soil conditions [40,41,42]. The soil CSC calculation method based on the Rhet data was also favored by a study that compared different soil CSC estimation methods based on chamber measurements against reference estimates using eddy covariance measurements. In this study, it was concluded that an approach based on Rhet data provided results that agreed better with the reference results, as compared to methods that use Rfloor measurement data [25]. It was found that, although both soil CSC estimation methods are sensitive to biases introduced by the soil C input and output data used in the calculations, the approach based on Rhet data was more applicable due to the relatively simple calculation approach of subtracting Rhet measurement results from the soil C input data, while complicated modeling of ecosystem photosynthesis and respiration is needed to calculate soil CSC using Rfloor measurement data.

To recalculate the study results of the Rfloor for individual study sites to Rhet, a factor ranging from 0.52 to 0.66 (mean 0.6) was used. Therefore, the calculated Rhet may be overestimated, as the equation applied to determine the Rhet/Rfloor recalculation factor was elaborated by comparing the data of Rhet and Rfloor, excluding aboveground autotrophic respiration (Rs), from studies conducted in both boreal and temperate zones [25]. However, such an assumption may be speculative as, according to the data compiled by more recent metanalytical reviews, the Rhet/Rs determined by the trenching method in boreal coniferous forests ranges from 0.36 to 1.03 (mean: 0.73, with a standard deviation of 0.18). Therefore, the approach used to calculate the Rhet in the study may be considered conservative from the GHG inventory perspective as underestimation of soil C loss is not favorable in the elaboration of soil emission factors.

The results of Rhet calculated for individual study sites with forest cover ranging from 2.9 to 4.4 t C ha⁻¹ year⁻¹ fall within the range of results of the Rhet estimated by direct measurement in other studies in boreal forests. The Rhet of forestry-drained peatlands reported in the results of a Finnish study ranges from 1.46 to 6.70 t C ha⁻¹ year⁻¹. [16]. The Rhet estimated in a 30-year-old Scots pine plantation (former cropland) with organic soil situated in the middle of a boreal climatic zone was 4.80 t C ha⁻¹ year⁻¹ [43], while the quantified Rhet of 12 afforested organic soil cropland sites and six cutaway peatlands in Finland ranged from 2.07 to 5.39 t C ha⁻¹ year⁻¹ and from 2.76 to 4.79 t C ha⁻¹ year⁻¹, respectively [44]. The results of another study carried out in Finland showed an average Rhet of 2.38 t C ha⁻¹ year⁻¹ in a pine-dominated drained mire. It was estimated that the annual Rhet of forestry-drained peatlands in central Estonia and southern and northern Finland ranges from 2.48 to 5.15 t C ha⁻¹ year⁻¹ [13]. The consistency of the Rhet estimated in our study with previously reported values indicates that the use of the Rfloor recalculation method is applicable for studies conducted in the hemiboreal zone.

4.2. Soil Respiration

In our study, the difference between the mean measured instantaneous Rfloor in drained sites (7.35 ± 0.89 t C ha⁻¹ year⁻¹) and undrained (7.02 ± 0.96 t C ha⁻¹ year⁻¹) forest stand sites were found to be insignificant (p = 0.34). However, the differences between the measured mean instantaneous Rfloor in forest stands (6.84 ± 0.56 t C ha⁻¹ year⁻¹) and clear-cuts (10.08 ± 1.96 t C ha⁻¹ year⁻¹) were significant (p = 0.002). The tendency of similar soil respiration in drained and undrained sites, as well as increased emissions in areas with no forest cover, was also observed in a previous study. For instance, in afforested lowland raised peat bogs in Scotland, it was found that Rfloor was slightly higher in drained sites (4.53 t C ha⁻¹ year⁻¹) compared to undrained sites (3.35 t C ha⁻¹ year⁻¹), while in undrained areas with no forest cover, the estimated Rfloor was 6.95 t C ha⁻¹ year⁻¹ [45].

In our case, the insignificant difference between average soil respiration in drained and undrained study sites is mainly determined by the absence of correlation between the measurement data of the Rfloor and groundwater level. The average groundwater level in the drained study sites (mean 55 ± 2 cm) was on average 18 ± 2 cm deeper than in the undrained sites (mean 35 ± 3 cm); however, the groundwater level measurement results have weak, insignificant correlation (r = 0.3) with soil respiration data. We were not able to find empirical reasons for having a significantly higher soil respiration rate in clear-cut study sites compared to sites with forest cover. By evaluating the data of Ts and atmospheric temperature measurements, we observed that there was no higher Ts sensitivity to variation in atmospheric temperature. Linear regression models of characterizing the relationship between Ts and soil temperature in forest stands and clear-cuts were not statistically significant. Therefore, we concluded that Ts response to changes in atmospheric temperature was not different in both study site groups and increased warming of soil in clear-cuts was not the reason for elevated soil respiration. Most likely, the increased emissions are induced by soil disturbances of mechanized harvesting [46] and the decomposition of harvesting residues [47].

The annual Rfloor in clear-cuts with peaty gley soil, estimated by the previous study as ranging from 6.5 ± 1.6 to 7.1 ± 1.7 t C ha⁻¹ year⁻¹ [48], which is similar to our estimation of the annual Rfloor in drained and undrained clear-cuts, i.e., 7.6 ± 2.2 and 7.9 ± 1.2 t C ha⁻¹ year⁻¹, respectively. Additionally, the estimated Rfloor in the study sites with forest cover, which ranged from 4.4 to 8.0 t C ha⁻¹ year⁻¹, is similar to the range of the Rfloor estimated in other studies (2.73 ± 0.55 to 5.18 ± 1.09 t C ha⁻¹ year⁻¹) conducted in boreal forests [49,50]. Furthermore, the mean Rfloor was found to be significantly higher in drained coniferous forests with organic soil (from 2.45 to 5.18 t C ha⁻¹ year⁻¹) than in undrained mire forest sites (from 2.18 to 3.27 t C ha⁻¹ year⁻¹), although the drained sites were all moist [51]. This may be in line with the observations made in our study that the groundwater level may have no significant impact on Rfloor. Additionally, in a study aimed at creating soil respiration prediction models, it was concluded that by adding the water table depth into the models as an explanatory variable, the goodness of fit was not improved and the prediction power was not statistically significantly improved [52]. Even though, in some cases, the average water table depth can be significantly correlated with annual respiration values in peatlands [53], soil temperature alone is generally sufficient to explain the variation in soil respiration. The reasons why groundwater level can be used as a Rfloor predictor only in some areas can be further studied.

4.3. Soil C Input by Litter

The study results of the average annual soil C input by LF in drained and undrained forests with nutrient-rich organic soils ranging from 1.6 ± 0.3 to 2.2 ± 0.4 C ha⁻¹ year⁻¹, respectively, are within the uncertainty range of the average observed values in the coniferous and deciduous forests of Northern Europe, 1.7 ± 1.1 and 1.5 ± 0.7 t C ha⁻¹ year⁻¹, respectively [54]. While similar relationship tendencies with the basal area have been recognized, higher estimated values of the average soil C input by LF of coniferous (1.82 ± 0.02 t C ha⁻¹ year⁻¹) and silver birch stands (2.07 ± 0.03 t C ha⁻¹ year⁻¹) with drained organic soils were found in a recent Latvian study [55]. This points out that the average soil C input values used in the estimates of forest C balance or comparison of litterfall biomass across different studies can lead to considerable inaccuracies. In our study, as well as in a previous local study [55], it is recognized that basal area provides the highest prediction power of litterfall biomass compared to other commonly used forest stand characteristics. Therefore, the variation in the LF data acquired in the study can be explained with the basal area of the forest stands studied. A limited number of study sites restricted the ability to compare the relationship between the basal area and LF in the drained and undrained sites separately.

The mean tree fine-root production, ranging from 0.6 ± 0.6 t ha⁻¹ year⁻¹ in drained deciduous forests to 1.0 ± 0.2 t C ha⁻¹ year⁻¹ in undrained coniferous forests with nutrient-rich organic soil, as estimated in this study, is significantly lower than those mostly reported in previous studies. The mean fine-root production reported for Northern Europe was 2.84 ± 1.52 t ha⁻¹ year⁻¹ in coniferous forests and 1.99 ± 1.01 t ha⁻¹ year⁻¹ in deciduous forests [56]. Lower estimated fine-root production values may be explained by methodological underestimation or the phenomenon whereby the growth of trees in nutrient-rich soil requires less biomass of the fine roots to ensure a sufficient intake of water and nutrients. Higher fine-root productivity in stands with less fertile soils has been noticed in several studies [22,57,58,59]; however, the opposite relationship has also been found [60]. Such assumptions may also be contradicted by the annual fine-root production in forests with drained nutrient-rich soil, which ranged from 1.81 to 3.02 t ha⁻¹ year⁻¹, as reported in an Estonian study [12]. Most likely, the reason for underestimation arises from methodology, as the average uncertainty of acquired results also ranges from 30 to 161 % (mean 71 %) in study sites with different dominant tree species and soil moisture regimes. The study period of one year was not sufficient for fine root production estimates by the ingrowth method, as one vegetation season was not sufficient to mitigate the effects of disturbance introduced by the installation of ingrowth bags.

The annual soil C input by GV has not been studied extensively, and the available results are often not comparable due to different study methods and the different vegetation components included in the calculations. In Estonia, the estimated primary production of dwarf shrubs and grasses reached 0.4 t C ha⁻¹ year⁻¹ in spruce stands and ranged from 0.6 to 1.0 t C ha⁻¹ year⁻¹ in pine stands [22]. By using the biomass of herbs and grasses, the prediction models elaborated by a study conducted in Finland [33], taking into account the age distribution of Latvian forests, the average weighted annual soil C input by aGV and bGV ranged from 0.34 ± 0.01 to 1.29 ± 0.202 t C ha⁻¹ year⁻¹ in birch and pine stands with drained organic soil, respectively [55]. The higher annual soil C input by herbs and grasses, which ranged from 0.6 to 3.2 t C ha⁻¹ year⁻¹ in forest stands and from 2.9 to 4.0 t C ha⁻¹ year⁻¹ in clear-cuts estimated in our study, can be explained by the forest site types characterized by nutrient-rich soils included in this study and the positive correlation found between GV biomass and soil fertility characteristics.

4.4. Annual Net Soil CSC

According to general opinion, the drainage of organic soil increases CO₂ emissions and reduces soil C stock; however, the results of previous studies on the effect of organic soil drainage on GHG emissions are ambiguous. The empirical data collected during this study shows that nutrient-rich organic soil in forest stands is a net CO₂ sink, but the soil in clear-cuts is a net CO₂ source. We estimated that during the study period, the average annual soil CSCs were 0.4 ± 0.4 t C ha⁻¹ year⁻¹ in undrained and 0.1 ± 0.4 t C ha⁻¹ year⁻¹ in drained forest sites, while in clear-cut estimated soil net C balance is −0.9 ± 0.7 and −0.4 ± t C ha⁻¹ year⁻¹ in drained and undrained sites, respectively. The observation of soil in forest stands acting as a C sink is in agreement with the conclusion reached in a previous local study on the CSC of drained moderately nutrient-rich forest soils [14] and can be explained by site productivity induced increased C input by litter that fully compensates soil C loss by respiration. The reason for soil in clear-cuts being a C source was increased C loss by respiration and reduced soil C input by litter. Although mean soil C input by ground biomass in clear-cuts (3.55 ± 0.37 t C ha⁻¹ year⁻¹) was considerably higher than in forest stands (1.65 ± 0.40 t C ha⁻¹ year⁻¹), that did not compensate for the increase in Rhet by average 0.8 t C ha⁻¹ year⁻¹ compared to forest sites and the absence of C input by litterfall (average 1.8 ± 0.5 t C ha⁻¹ year⁻¹) and the fine roots of trees (0.4 ± 0.2 t C ha⁻¹ year⁻¹).

In addition to calculating the annual soil CSCs as a sum of soil C balance components quantified in a monitoring year of the study, we modeled an inter-annual soil CSC within a 240-year forest management period by using variables of stand characteristics as predictors. As a result, by taking into account the impact of forest stand development (age and basal area) and the stages of forest land (forest stand and clear-cut), the average soil CSC with high variability was estimated. The range of estimated annual soil CSCs is mainly determined by two soil C balance components with high inter-annual variability determined by the stage of forest stand development—LF and ground GV. The most significant impact on variability is introduced by GV, which determines on average 63% of the total soil C input by litter, while LF impacts on average 24% of annual soil C input. Empirical data with the highest uncertainty is soil C input by FR, aGV, and bGV with confidence intervals of 71%, 41%, and 37%. The annual soil C input by GV ranges from an average of 1.68 ± 0.36 to 4.88 ± 1.40 t C ha⁻¹ year in clear-cuts and forest stands with various dominant tree species signifying the importance of reducing the uncertainty of annual soil CSC estimations introduced by data on GV.

The results of the study are in line with the previous studies, which showed that the soil C stock does not change and can even increase after the drainage of organic soil in boreal forests [11,24,33,34]. It can be expected that soil respiration may be considerably higher in nutrient-rich site types compared to site types with less fertile soils [13], which may lead forests with fertile drained organic soils to be a source of CO₂ emissions for the following reason. The Rfloor of forestry-drained peatlands estimated in Finland showed a clear diminishing trend in annual soil respiration from the most to the least fertile site types, and ranged from 3.8 to 12.10 t C ha⁻¹ year⁻¹ [16]. However, the results of this study indicate that both undrained and drained nutrient-rich organic soil in forest stands can still be a net C sink. The differences in the calculated annual soil CSC across various studies may be due to the different methods applied and the inclusion of different soil C input components in the calculations [52], as well as the uncertainty of these values, since the data available, especially for belowground litter, are highly uncertain, most often due to difficulties in acquiring such data [56]. For example, drained forest peatlands were identified as a CO₂ source also in Sweden; however, the estimated annual soil C loss of −2.29 t C ha⁻¹ year⁻¹, calculated by subtracting Rhet from the soil C input [61], is considerably higher than that in our study, whereas in Finland, peatlands drained for forestry were found to be a net CO₂ sink (removals from 0.2 to 0.252 t C ha⁻¹ year⁻¹) estimated by the soil C inventory method [62]. In our study, the inclusion of forest land status as a clear-cut estimation of annual soil CSC determined if nutrient-rich drained organic forest soil acts as a CO₂ sink or source. The differences highlight the importance of harmonizing soil CSC estimation methods to improve the comparability of country-level GHG inventory results.

The results of soil CSC acquired in this study can be further improved by both more extensive studies and by conducting direct Rhet measurements or evaluating the proportion of Rhet/Rfloor under national conditions. Instead of using static annual soil C input value, approaches to model inter-annual litter biomass variations based on forest stand variables and climatic conditions should be elaborated for the estimation of soil C balance by offsetting the annual Rhet modeled using the annual data of air or soil temperature. The inter-annual variation in hourly or diurnal temperature data may have a significant impact on modeling soil respiration using previously elaborated equations characterizing the relationship between soil respiration and soil or atmosphere temperature. The choice of using temperature data of one year or time period characterizing climate, as well as the use of daily mean or hourly mean temperature data, may have a considerable impact on soil respiration modeling results, which should be considered in future studies.

5. Conclusions

The drained and undrained nutrient-rich organic soils in the forest stands monitored for one year in this study were a CO₂ sink, while the soil in clear-cuts acted as a CO₂ source. The soil in clear-cuts acting as a CO₂ source was determined by increased soil respiration rates and the absence of soil C input by litterfall and the fine roots of trees. The significantly increased soil C input by ground vegetation in clear-cuts mitigated this effect; however, the significantly increased soil respiration rate and reduced soil C input by other sources were not fully compensated. If forest management cycles are considered, including forest land state as a clear-cut, drained nutrient-rich organic soil in managed forests is a CO₂ source, while the soil C stock increases in undrained soil, according to the methodology applied in the CSC calculations.