The impact of wildfire on biogeochemical fluxes and water quality in boreal catchments

Wildfires are the major disturbance in boreal ecosystems and are of great importance for the biogeochemical cycles of carbon (C) and nutrients. However, these fire-induced impacts are hard to quantify and are rarely assessed together at an ecosystem level incorporating both aquatic and terrestrial environments. Following a wildfire in Sweden in an area with ongoing monitoring, we conducted a pre-fire (9 years) and post-fire (4 years) multi-catchment investigation of element losses (combustion and leaching) and impacts on water quality. Direct C and nitrogen (N) losses through combustion were ca. 4500 and 100 g m⁻², respectively. Net CO₂ loss associated with soil and biomass respiration was ∼ 150 g C m⁻² during the first year, but the ecosystem started to show net CO₂ uptake in June 3 years post-fire. Aquatic C and N losses the first 12 months post-fire were 7 and 0.6 g m⁻², respectively. Hence, soil respiration comprised a non-negligible part of the post-fire C loss, whereas aquatic C losses were minor and did not increase post-fire. However, other elements (e.g. Ca, S) exhibited ecologically relevant increases in fluvial export and concentration with large peaks in the immediate post-fire period. The temporal dynamics of stream concentrations (Ca²⁺, Mg²⁺, K⁺ ,SO${}_{\mathrm{4}}^{-\mathrm{2}}$, Cl⁻ ,NH${}_{\mathrm{4}}^{+}$, total organic N) suggest the presence of faster- and slower-release nutrient pools with half-lives of around 2 weeks and 4 months which we attribute to physicochemically and biologically mediated mobilization processes, respectively. Three years after the fire, it appears that dissolved fluxes of nutrients have largely returned to pre-fire conditions, but there is still net release of CO₂.

1 Introduction

Wildfires are the major disturbance agent in boreal ecosystems and are expected to increase in size and frequency (Flannigan et al., 2009). Wildfires have a large impact on biogeochemical cycles, and emissions of CO₂ to the atmosphere from more frequent and larger wildfires could generate a positive climate feedback unless the carbon (C) emitted is swiftly re-sequestered (Bond-Lamberty et al., 2007; Smithwick et al., 2005). Wildfires also influence the biogeochemical cycles of nitrogen (N) and major cations (Brais et al., 2000; Grier, 1975; Smithwick et al., 2005), which can influence post-fire ecosystem productivity, an issue which has been discussed for decades (e.g. Ahlgren and Ahlgren, 1960; Grier, 1975). Losses occur both as emissions during the fire and through post-fire losses via runoff. However, these fire-induced impacts are hard to quantify and are rarely assessed at an ecosystem level including both aquatic and terrestrial environments (Amiro et al., 2010; Brais et al., 2000; Rhoades et al., 2019; Turner et al., 2007). Comparing post-fire responses to pre-fire conditions is also problematic because wildfires rarely take place at locations with pre-fire measurements. Here we present a unique pre- and post-fire multi-catchment investigation of water quality and element cycling in boreal Sweden.

Boreal wildfires often consume a large portion of the fuel in the form of ground vegetation and can also consume the upper organic soil (Amiro et al., 2000; Turetsky et al., 2011). Up to 90 % of the emitted carbon typically comes from the organic soil layer, and in North America, such C emissions are estimated to be on average 3000–4000 g C m⁻² (Turetsky et al., 2011; Walker et al., 2018). In drained peatlands, the increased exposure of organic soil to oxygen means that C losses can be 1 order of magnitude larger than uplands and undrained peatlands (Granath et al., 2016). In addition to C, N is also emitted in large quantities during fires (Johnson et al., 2007) as it starts to volatilize at 200 ^∘C (Knicker, 2007). This contrasts to other nutrients (e.g. K, P) that require a combustion temperature above 760 ^∘C (Knicker, 2007), which rarely occurs. Although N losses can potentially influence long-term ecosystem productivity (Tamm, 1991), few studies have quantified N emissions via this pathway (Brais et al., 2000; Johnson et al., 2007). Studies that have quantified ecosystem C and N emitted during wildfires are still scarce and are lacking for northern Europe, impeding our understanding of how wildfires alter major geochemical cycles.

Boreal wildfires do not only cause direct emissions of C and nutrients but can also alter their fluvial transport and thus downstream water quality (Bladon et al., 2014). To what extent this is true for C does, however, depend on the compound measured, catchment characteristics, and probably fire severity (Santos et al., 2019). Studies have shown negative, little, or no effect on the total amount of dissolved organic carbon (DOC) exported post-fire (see discussion in Evans et al., 2017; Rodríguez-Cardona et al., 2020), whereas DOC aromaticity and particulate organic carbon (POC) export can increase (Burd et al., 2018; Evans et al., 2017; Olefeldt et al., 2013). More striking is the increase in available macronutrients and other elements that are released from the burned organic top layer. Typically, the loss of soil cation exchange capacity resulting from the combustion of organic soil, together with the combustion of biomass, leads to the release of exchangeable cations (e.g. Ca²⁺, Mg²⁺, and K⁺; González-Pérez et al., 2004). These ions are easily exported to streams and lakes and can lead to an increase in runoff pH. On the other hand, many studies have shown post-fire peaks in sulfate (SO${}_{\mathrm{4}}^{\mathrm{2}-}$), chloride (Cl⁻), and nitrate (NO${}_{\mathrm{3}}^{-}$) due to a combination of release from soil and reduced biological demand (notably for NO${}_{\mathrm{3}}^{-}$) (Bayley et al., 1992; Bladon et al., 2008; Carignan et al., 2000; Lydersen et al., 2014; Mast and Clow, 2008). If acid anions (NO${}_{\mathrm{3}}^{-}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$ and Cl⁻) dominate over base cations, an acidity effect is observed in downstream waters (Lydersen et al., 2014). This acidification effect is enhanced in areas which have higher concentrations of stored S or N from historic deposition or have a high proportion of peatlands (Bayley et al., 1992; Evans et al., 2017). Lower pH increases dissolved P in the post-fire soil (Certini, 2005) and a long-term (3–5 years) increase in exported P in burned catchments has been reported across boreal Canada (Burd et al., 2018; Burke et al., 2005; Lamontagne et al., 2000; Silins et al., 2014). However, a high base cation concentration may counterbalance the downstream acidity effect (Carignan et al., 2000).

Nitrogen levels in runoff water normally increase dramatically post-fire (e.g. Bladon et al., 2008; Carignan et al., 2000). Following fire, soil organic nitrogen is either volatilized or converted into ammonium (NH${}_{\mathrm{4}}^{+}$), while nitrate (NO${}_{\mathrm{3}}^{-}$) is mainly formed from NH${}_{\mathrm{4}}^{+}$ through nitrification, a process which can continue for several years after the fire (Certini, 2005). With the loss of vegetation after a severe fire and limited potential for microbial immobilization due to a shortage of labile carbon, ammonium and nitrate cannot be retained within the ecosystem and are commonly leached out (Smith et al., 2011). Nitrate concentrations peak shortly after the fire, but the return time to reference values seems to vary from 2 to 9 years post-fire (e.g. Bladon et al., 2008; Carignan et al., 2000; Evans et al., 2017; Hauer and Spencer, 1998; Mast and Clow, 2008). In contrast to NO${}_{\mathrm{3}}^{-}$, NH${}_{\mathrm{4}}^{+}$ is expected to be held by the soil to a higher degree because it adsorbed onto negatively charged surfaces of soil particles (Mroz et al., 1980). However, a study observed NH${}_{\mathrm{4}}^{+}$ pulses that lasted over 2 growing seasons (Grogan et al., 2000).

Variation in surface water quality and fluvial transport in a boreal catchment is mainly controlled by landscape heterogeneity (Humborg et al., 2004). For example, the proportion of peatlands in a catchment has a major influence on surface water DOC and NO${}_{\mathrm{3}}^{-}$ which affect runoff pH through the release of organic acids (Buffam et al., 2007; Sponseller et al., 2014). Peatlands naturally retain sulfur under waterlogged conditions (in reduced organic forms and sulfides), so wildfires may lead to particularly high SO${}_{\mathrm{4}}^{\mathrm{2}-}$ leaching when peatlands burn. Beside peatlands, lakes upstream can act as buffers in the system by increasing residence time. This will dampen the water quality response to wildfire at the catchment outlet and possibly reduce the biogeochemical signal via element retention (e.g. in sediments). Despite the clear effect of landscape characteristics on water chemistry, we currently know little about what determines the magnitude or temporal dynamics of post-fire element leaching at the landscape scale.

In 2014, a large wildfire affecting established monitoring sites in Sweden created the opportunity to study ecosystem-level effects of wildfire on biogeochemical cycles in a managed boreal landscape. Whole-catchment studies are important in ecosystem science (Likens et al., 1970) but difficult to conduct at a detailed level, particularly in relation to unpredictable events such as wildfires. In our study, the burned area (circa 13 000 ha) consists of multiple catchments, allowing us to investigate local variation in post-fire responses. One of the catchment streams and one lake are included in the Swedish national water monitoring network, enabling comparisons with pre-fire data and with longer-term trends in water chemistry. Hence, compared to most studies, our study does not rely on a single catchment or only post-fire data (see Betts and Jones, 2009; Evans et al., 2017; Mast et al., 2016, for other before and after studies). In addition, it is rarely possible to study biogeochemical processes during the critical period immediately following a fire due to limited access to the area, as well as resource constraints.

The overarching aim of this study was to examine the impact of wildfire on element fluxes and water quality in boreal forests. Our first objective was to determine C and N losses through combustion during the fire and investigate how important these losses are compared to pre-fire soil pools, post-fire hydrologically exported C and N, and post-fire terrestrial C balance and plant regrowth. Secondly, we tested if hydrologically exported amounts of S, Ca, and K increased over the first 3 years post-fire and if such losses can affect their long-term soil pools. Our third objective was to quantify the magnitude and shape of the early post-fire flush and multi-year trend of nutrients in five streams and one lake. Here we also determined the form of element concentration decay curves (single or double exponential decay curves; Minderman, 1968) to understand post-fire biogeochemical cycling and ecosystem recovery.

2 Materials and methods

2.1 Study area

The study area is boreal forest located in southern Sweden (59^∘54^′50^′′ N, 16^∘09^′50^′′ E). It is located about 75 to 150 m above the sea level that has a low relief but is topographically complex. Between 1987 and 2016 the mean annual temperature was 6 ^∘C (January −3.3, July 17 ^∘C), and the annual precipitation was 687 mm. The forest is intensively managed using clear-cutting, planting, and thinning operations that create a mix of even-aged forest stands from recently cut areas to mature stands (> 100 years). Tree cover is dominated by Pinus sylvestris (particularly the catchments investigated here), shrub layer by Vaccinium myrtillus, V. vitis-idea, Calluna vulgaris, and Rhododendron tomentosum, and ground layer by Pleurozium schreberi, Hylocomium splendens, Polytrichum sp., and Cladonia sp. (see Gustafsson et al., 2019 for more details about the area). The area contains many small lakes (residence times mostly between 1 and 3 months) and has a high peatland coverage (10 %–35 %; Table 1, Fig. 1). The mineral soil consists of granitoid till and is general thin where peatlands are not present. A wildfire started on 31 July 2014 and burned over 12 d covering an area of ca. 13 000 ha. The fire was low intensity during the first days but spread rapidly when the wind speed increased and changed direction, and it became a high-intensity stand replacing fire across all catchments investigated in our study. Due to the high intensity, fire fighting efforts were mostly restricted to protecting populated areas. Half of the burned area was salvaged logged during the first year after the fire, while the other half was protected and left for natural regeneration.

Figure 1 Map over the burned area showing the sampled catchments, sampling points (terrestrial C and N loss and reference plots, as well as water sampling stations), and placement of the eddy-covariance towers. Most of the burned, non-peatland areas of the catchments Ladängsbäcken, Myckelmossbäcken, and Märrsjöbäcken/Märrsjön were salvage logged. Table 1 contains information about the catchments.

The burned area consists of multiple catchments. We defined five major catchments in ArcGIS 10.3 (ESRI, Redlands, USA ) by using the Swedish elevation model (resolution 2 × 2 m and elevation accuracy of 0.5 m; Lantmäteriet, 2014). When rain hits the surface it will run in the steepest slope direction which is determined in the elevation model. We delineated watersheds by grouping the surfaces of the steepest slopes with the same direction. Two of these catchments are within the perimeter of the nature reserve with little salvage logging (Gärsjöbäcken and Vallsjöbäcken), while two are largely salvaged logged (Myckelmossen and Märrsjön; Table 1, Fig. 1). All catchments were close to completely burned, and their outlets were placed just outside the burned area where water sampling were performed.

2.2 Stream water sampling and chemical analyses

Over 3 years post-fire we sampled outlet stream water from the five catchment outlets and near surface water from one lake (Märrsjön). One stream (Gärsjöbäcken) and the lake (Märrsjön) are included in the Swedish long-term monitoring programme (Fölster et al., 2014) and therefore have a long period of pre-fire data (something which is relatively rare in studies of wildfire impacts). We extracted data for the sites from 2005 to the present day (Miljödata-MVM, 2019). For all sites, post-fire stream sampling begun 2–3 weeks after the fire (ca. 1 week after the first major post-fire rain event, > 20 mm) and continued with high temporal resolution during the first 4 months and thereafter with longer intervals depending on season and stream. The lake was sampled slightly less frequently. The water sampling and subsequent water chemistry analysis were made according to the Swedish monitoring programme using standard methods at the SWEDAC-accredited (Swedish Board for Accreditation and Conformity Assessment) geochemical laboratory at the Department of Aquatic Sciences and Assessment at the Swedish University of Agricultural Sciences. Metal ions were analysed with inductively coupled plasma mass spectrometry (ICP-MS), and SO₄ and Cl were analysed by ion chromatography. NH${}_{\mathrm{4}}^{+}$ and NO${}_{\mathrm{2}}^{-}+$ NO${}_{\mathrm{3}}^{-}$ were analysed with an auto-analyser. Total organic carbon (TOC) and total N (TN) were analysed by combustion on unfiltered water samples (Shimadzu TOC-VCPH with a TNM-1 module). By using unfiltered water samples we include organic material that was washed out by erosion. In these boreal ecosystems the composition of TOC is completely dominated by DOC (Laudon et al., 2004). DOC was measured, together with TOC, in one stream during the first year, and these variables were highly correlated (r=0.98). We therefore use TOC as a proxy for DOC. Total organic N (TON) was calculated as follows: TON = TN − (NH${}_{\mathrm{4}}^{+}$)-N − (NO${}_{\mathrm{2}}^{-}+$ NO${}_{\mathrm{3}}^{-}$)-N.

2.3 Pre-fire soil conditions and carbon and nitrogen losses

We estimated shrub, moss, and organic soil C and N losses in the two largest catchments (Vallsjöbäcken and Gärsjöbäcken). Our large-scale sampling was based on a systematic 300×300 m grid. At each intersection of the grid, a 314 m² circular plot (r=10 m) was established for sampling (i.e. 300 m between each sampling plot). Within the plot we established two perpendicular transects with 41 sampling positions (every metre and in the centre). The high sampling density was chosen as burn severity is known to be extremely heterogeneous and spatial autocorrelation of organic soil depth is likely somewhere between 0.85 and 2.85 m (Kristensen et al., 2015). At each position, we registered the fire effect on the shrub layer (intact, only singed, only charcoaled stumps remaining, or totally consumed). For non-peaty soils (< 30 cm of organic matter), we measured the depth of the remaining soil organic layer (to nearest half centimetre) and recorded whether the top layer (moss–lichen + O_i horizon) had been consumed or not at each of the 41 positions within the plot. The ash layer (defined as “the particulate residue remaining, or deposited on the ground, from the burning of wildland fuels and consisting of mineral materials and charred organic components”; Bodí et al., 2014) was considered as remaining soil and was generally thin (0–0.5 cm). By including the ash layer in our measurements of remaining organic soil, we introduce additional uncertainty to our carbon loss estimates if C density is much different in this layer. To evaluate this effect we performed sensitivity analyses using ash C content, thickness, and weight from another study from the same burned area (Perez-Izquierdo et al., 2020). The plot mean was used to estimate depth of burn (DOB) as the predicted organic soil layer depth (based on reference sampling outside the burned area) minus the remaining depth (e.g. Kelly et al., 2016; Turetsky et al., 2011). In peatlands, we measured DOB at each position by measuring the distance between the post- and pre-fire positioning of the organic layer. We reconstructed the pre-fire position using the positioning of adventive roots on the basal area of tree trunks, positioning of horizontal tree roots, and positioning of remnants of the ground vegetation and peat mosses (for a detailed description of the methods see Kelly et al., 2016; Turetsky et al., 2011).

Table 1 Overview of the burned catchments, their land characteristics, and annual mean outflow water age (2014 August–2015 July). The last catchment (Märrsjön) is a lake catchment. Proportion logged is based on estimated salvage-logged area during the first year after the fire. Longitude and latitude (WGS84) indicate the sampling location. See also the map in Fig. 1. Long-term monitoring catchments are indicated with an asterisk (*).

Carbon and nutrient losses during the fire were estimated for the organic soil layer and ground vegetation. To do this we needed to reconstruct the pre-fire organic soil thickness, bulk density, and nutrient content (C, N, S, K, Ca, P) of the organic soil layer, moss and/or lichen layer, and ground-layer cover of shrubs to calculate their biomass and ultimately their C and N content. Using the same protocol as for the burned plots, we collected data from 10 reference transects in the unburned surroundings, amounting to up to 57 plots (Fig. 1). These transects were placed from hilltops to valley with five to seven plots per transect, covering young to old forests, similar to the area burned. Peatlands were not included as we estimated depth of burn directly in these habitats.

For reference data on the organic soil layer, we sampled three to five soil cores (d=10 cm, depth = 5–30 cm depending on terrain) per plot and split them into a living moss and/or lichen section including the O_i horizon and a decomposed section (O horizon consisting of horizons O_e and O_a). Each section was dried (65 ^∘C, until no further weight loss occurred), weighed, mixed, and thereafter analysed for total element mass by Forest Research, UK. Elements were measured on a mass basis (g kg⁻¹) and converted into element bulk density (BD; g cm⁻³). We used the DOB estimates and bulk density values (moss–lichen layer + O_i and O_e+a horizon) to calculate the soil C and N losses per area (DOB × BD). Unburned reference sites have often been used as controls to estimate fire-generated C and N losses (e.g. Kelly et al., 2016; Turetsky et al., 2011) and produce estimates similar to studies that used both pre- and post-fire measurements (Johnson et al., 2007). For peatlands we used published data on BD (5 cm depth interval; Granath et al., 2016) for boreal drained and undrained peatlands as the treed peatlands in the burned area in general are drained. Peat C content and N content were assumed to be 55 % and 2 %, respectively (Minkkinen and Laine, 1998).

In our study we call these losses for direct losses (or emissions), meaning that they were predominantly lost from the soil and ground vegetation at the time of the fire. DOB data were collected within 1 year post-fire, and for uplands they were based in the remaining organic soil layer. Hence, there is a possibility that we include other early losses (e.g. fluvial and respiration losses) in our upland direct emission estimates.

We estimated ground vegetation cover in the reference plots by recording the presence/absence of dwarf shrubs at 41 positions within each plot. To convert cover to biomass we used species-specific relationships between cover and biomass for the major shrubs species (Vaccinium myrtillus, V. vitis-idaea, Calluna vulgaris, and Rhododendron tomentosum). In a second step, we scaled up C and N losses to catchment level by using the average losses for upland and peatland weighted by their coverage, respectively. Peatland cover was retrieved from the Swedish Geological Survey database (https://apps.sgu.se/kartvisare/kartvisare-torv.html, last access: 1 February 2020).

C and N losses from standing trees were not estimated. It is very hard to make reliable quantifications of such losses (amount of fine branches and needles consumed), and the fuel amount varies with stand density and age. A typical pine stand in the burned area may have 750 stems per hectare and a stem diameter between 15 and 20 cm and be 15–20 m high. This gives about 0.5 kg m⁻² C stored in living branches and needles and 0.15 kg m⁻² C only in needles (calculated using allometric equations from Marklund, 1988). Only 21 % of the area experienced 100 % crown damage and about 50 % between 50 % and 100 % damage (Gustafsson et al., 2019). Charred needles and fine branches were still visible in the burned pine crowns, indicating small losses from the trees and likely amounting up to a few per cent of the total C loss in forested areas.

2.4 Measuring CO₂ fluxes

Net ecosystem exchange (NEE) of CO₂ was measured by eddy covariance (EC) at two locations within the burned area (Fig. 1). Each EC system comprised a CSAT3 sonic anemometer and an EC155 closed-path gas analyser as an integrated system (CPEC200, Campbell Scientific, Logan, UT, USA). The sensors were mounted on a boom at the top of a 2 m tripod. Measurements were made at 10 Hz using a CR3000 datalogger (Campbell Scientific, Logan, UT, USA). Meteorological measurements including air temperature, solar radiation, and soil moisture and temperature at 5 cm depth were recorded at the same location as 30 min averages. Raw 10 Hz EC data were aggregated to calculate 30 min average CO₂ fluxes, and overall fluxes were calculated according to the EUROFLUX methodology for error correction and gap-filling (Aubinet et al., 1999; Lee et al., 2004). In particular, detrending was applied using a digital recursive filter with a time constant of 2000s, and the covariance matrix was aligned with the mean wind vector by a two-fold coordinate rotation on a half-hourly basis. Data analysis was done using R (R Development Core Team, 2016) and the R package openair (Carslaw and Ropkins, 2012). The EC systems were installed in April 2015 due to limitations in accessing the burned area, and CO₂ fluxes prior to that date (autumn–winter) were modelled. For a more detailed description of the data processing and gap-filling techniques used, see Hadden and Grelle (2017).

2.5 Element budget calculations

To make approximate element budgets we combined estimates of pools and fluxes in the system. Our aim was not to make a complete budget but rather to contrast immediate changes in stocks (assumed to be direct gaseous emissions for N and C) during the fire and subsequent (leached out or net ecosystem CO₂ exchange) losses from the ecosystem. Pre-fire element pools were derived from reference sites, and emissions were estimated from DOB (see text above). This was done for the two major catchments (Gärsjöbäcken, Vallsjöbäcken) for which we had DOB measurements. Fluvially transported material was calculated based on stream flow and water element concentrations. Flow data were based on S-HYPE (Strömqvist et al., 2012), the national application of the HYPE hydrological model (Lindström et al., 2010). HYPE is a process-based daily time-stepping catchment model. In a HYPE model application the modelled domain is divided into sub-basins with unique distributions of hydrological response units (HRUs). These HRUs are typically a combination of specific land uses and soil types. The soil profile of each HRU may contain up to three soil layers. Runoff of water from the soil layers including overland flow are simulated and summed for each HRU and routed through the network of rivers and lakes in the model. The Vallsjöbäcken catchment was extracted from the national model application and calibrated against local pre-fire and post-fire streamflow data using an automatic calibration routine. Pre-fire data were obtained from a stationary streamflow gauging station in operation until the early 2000s. Post-fire streamflow time series were derived from data from installed pressure transducers and a rating curve developed from the recorded water level and flow measurements. The post-fire model was validated against streamflow data derived from the transducer installed in Gärsjöbäcken. Using this model we also extracted daily estimates of the average residence time of water in the drainage network upstream of the sampling point.

Element mass flow was calculated as daily flow times element concentration. As element concentration was not measured daily we used predicted values from a model that made linear predictions between time points. This approach (period-weighted) was chosen over a model based on flow–concentration relationships because such relationships were weak in our data, indicating that non-hydrological factors dominated observed temporal variations (see Results). Our approach is recommended by Aulenbach et al. (2016) when there is a weak concentration – discharge relationship and the load estimate error should not be larger than 5 %–10 % (Aulenbach et al., 2016). Element outflow was aggregated over time, and we present values for 3 years pre-fire (for Gärsjöbäcken catchment, the long-term monitoring site) and for 3 years post-fire (Gärsjöbäcken and Vallsjöbäcken).

2.6 Element decay curves and pH modelling

For solutes that showed a single “pulse” response to the fire (Cl⁻, Ca²⁺, Mg²⁺,K⁺ ,SO${}_{\mathrm{4}}^{-\mathrm{2}}$, NH${}_{\mathrm{4}}^{+}$, TN), we fitted exponential decay curves to observed concentrations in order to derive a set of diagnostic parameters describing the magnitude of fire response and rate of recovery to pre-fire baseline conditions. This procedure was undertaken at the four streams with sufficient data to support curve fitting: Myckelmossbäcken, Ladängsbäcken, Gärsjöbäcken, and Vallsjöbäcken. Based on an initial assessment of the data, it was apparent that some solutes did not follow a simple (single) exponential decay curve, whilst in all cases solute concentrations converged on a non-zero baseline concentration towards the end of the measurement period. Therefore we conceptualized the change in solute concentrations according to Eq. (1):

$$\begin{array}{}\text{(1)}& \begin{array}{rl}{C}_t& =C_{\mathrm{baseline}}+C_{\mathrm{fast}}\times {\mathrm{0.5}}^{\left(t/t\mathrm{1}/{\mathrm{2}}_{\mathrm{fast}}\right)}\\ & +C_{\mathrm{slow}}\times {\mathrm{0.5}}^{\left(t/t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}\right)},\end{array}\end{array}$$

where C_t represents solute concentration at time t, C_baseline is the average concentration of a solute in the absence of fire effects, and C_fast and C_slow are the maximum post-fire concentrations of two exponentially declining pools with associated half-lives of $t\mathrm{1}/{\mathrm{2}}_{\mathrm{fast}}$ and $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$, respectively.

For each solute at each site, we fitted non-linear decay curves (Eq. 1). First, we located the time of peak measured concentration at each site (which was not necessarily the same at all sites nor was it the first measurement post-fire) as time zero. Next, we estimated C_fast, C_slow, $t\mathrm{1}/{\mathrm{2}}_{\mathrm{fast}}$, and $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$ for each solute time series by using a Bayesian approach in the R package brms version 2.10 (Bürkner, 2017). To regularize estimation we used weakly informative (proper) priors based on expected values: mean and SD 10 for C_baseline and $t\mathrm{1}/{\mathrm{2}}_{\mathrm{fast}}$ and mean 100 and SD 25 for C_fast and C_slow and $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$. A least-square estimation gave similar results but was sensitive to starting values for each model. The fast pool was tested by examining if the 95 % credible intervals of the fast-pool parameters included zero. With the fitted models we defined pool half-lives; the amount and relative proportion of peak measured concentrations associated with baseline and fast- and slow-decay pools; and the ratio of peak to baseline concentrations for each site and solute combination.

We modelled pH and charge of organic anions (RCOO⁻) following the approach by Köhler (2000) which is based on TOC, alkalinity, and pCO₂=2.8 using the CBALK approach. Charge balance with respect to buffering capacity and organic anions is achieved through iteration until a charge balance criterion of positive and negative charges (< 0.1 µeq L⁻¹) is met. The pH measurements were taken coincident with the water samples to validate this model.

2.7 Leaf area index

To examine post-fire plant regrowth, we extracted remotely sensed leaf area index (LAI) at peak growing season (15 June–28 July) for 2014 (before fire) to 2019. We downloaded MODIS LAI data (product: MCD15A2H) with a 500 m pixel size and 8 d averages (Myneni et al., 2015). We filtered out “bad” pixels using the quality layers (e.g. pixels with clouds and high aerosol content). Pixels covering more than 25 % water were also removed from further calculations. Finally, we extracted the mean values for each catchment and year. MODIS data were downloaded in R using the MODISTools package (version 1.1.1, Tuck et al., 2014), and calculations were performed with the raster package (version 3.0-7; Hijmans et al., 2019).

3 Results

3.1 Element losses and C fluxes

C and N losses from the soil and ground vegetation during the fire (assumed to be emissions) were similar in the two focus catchments (Table 2). In forest (non-peaty) soils and ground vegetation, most of the C and N losses were from the O horizon, while the contribution of the shrub vegetation was negligible (ca. 2 %). The moss and/or lichen layer (pre-fire thickness 28 ± 15 mm, mean ± SE) was, with rare exceptions, completely consumed by the fire. On average, 12 mm of organic soil remained after the fire (compared to an estimated 98 ± 53 mm pre-fire), and the organic soil C and N stock had been drastically reduced (−88 %).

Table 2 Emissions during the fire and hydrologically exported masses for major elements. Pre-fire data are estimated for one catchment (Gärsjöbäcken catchment) as a 4-year average (2010–2013). Post-fire exported masses are calculated over 12 months periods after the fire (August–July) for 3 years (2014–2017) with an expected error of 5 %–10 % (Aulenbach et al., 2016). Amounts are calculated from total N, total P, SO${}_{\mathrm{4}}^{\mathrm{2}-}$-S, Ca, K, and total organic C (TOC in streams). Emission is estimated losses to the air during the fire, excluding trees (dead and living). Pre- and post-fire storage is only for the forest habitat (excluding peatlands). Water is discharge per year, normalized for catchment area, and expressed in metres.

^a Average loss over peatlands and forests; ^b peatlands excluded; ^c for 2014, 2015, and 2016.

Fluvial element transport was controlled mainly by element concentration as we found no evidence that element concentration was a function of stream flow. For the two catchments, flow explained at the most (for K at Vallsjöbäcken) 17 % of the variation in element concentration, followed by SO${}_{\mathrm{4}}^{\mathrm{2}-}$ with 10 %–11 % explained variation for the two catchments (Fig. S1). In the Gärsjöbäcken catchment that had pre-fire data, the streamflow and element concentration relationship was equally weak the years before the fire (R² < 20 %). Furthermore, a pre- vs. post-fire comparison showed that fluvial losses increased drastically for all elements and were around 5 times higher during the first year, except for S that was 26 times higher (Table 2). In the third year post-fire, S and P still showed higher values than before the fire, whilst Ca and K had returned to pre-fire levels. For Vallsjöbäcken catchment, fluvial losses were overall lower than for Gärsjöbäcken, but the temporal trend was almost identical. Discharge was substantially higher the first year (50 %–60 %) in the two catchments but thereafter similar to the pre-fire values.

Carbon fluxes were similar at the two sites, and, on average, these two sites lost 158 g C m⁻² the first year and in total ∼ 440 g m⁻² (426 and 456 g m⁻²) over 3 years (Fig. 2a, b). This is about 10 % of the C lost in the fire. Merging all C losses and fluxes over the first 3 years, we estimated the total C loss to be circa 4900 g m⁻² in the two catchments. There was a net C loss for all months except for a few summer months close to 3 years post-fire. This trend towards a net carbon uptake was mirrored in the large-scale vegetation regrowth data. Regrowth (here as LAI) occurred at a similar rate among the burned areas of the catchments (Fig. 2c). Ladängsbäcken, where 28 % of the catchment area did not burn, showed a weaker response when LAI was estimated for the whole catchment (lowest value 1.84).

Figure 2 Carbon fluxes from two burned areas measured by eddy covariance and changes in leaf area index (LAI) for the five catchments. (a) Monthly CO₂-C fluxes and (b) cumulative CO₂-C flux. The dashed line indicates the period of modelled fluxes. Positive values mean emissions to the atmosphere. (c) Changes in mean summer (June 15–July 28) LAI over time of the burned parts of the catchments. First year (2014) shows LAI prior to the fire. Data are from MODIS (500 m grid). The fire occurred in August 2014 and is indicated by a dashed vertical line.

3.2 Water quality

Nitrate and ammonium concentrations increased rapidly post-fire, and ammonium quickly decreased and stabilized within 12 months in all catchments (Fig. 3). Nitrate, however, continued to show spring pulses. Soluble P also increased in streams, but the magnitude varied, and there are indications of winter–spring pulses. SO${}_{\mathrm{4}}^{\mathrm{2}-}$, Ca²⁺, and K⁺ concentrations followed the same pattern as ammonium and had stabilized after a year, except for K that returned at a slower pace. Fire had a marginal effect on pH and TOC in streams (Fig. 3; Supplement Figs. S2, S3). A short acidification pulse (0.5–1 pH unit) occurred during the first few months, but then pH slowly increased over time. Analyses of ions indicate that the pH was relatively stable after the fire because increases in acidity caused by SO${}_{\mathrm{4}}^{\mathrm{2}-}$ were counterbalanced by organic acids and an increase in base cations (Ca²⁺, Mg²⁺, Na⁺, and K⁺) (Supplement Figs. S2–S4). The pH modelling exercise resulted in a median difference between measured and modelled pH in this data set of 0.19 pH unit. The large majority (> 90 %) of the measured pH could be modelled within 0.5 pH units, which is in line with earlier similar studies (Fig. S5).

Figure 3 Temporal changes in concentration of nutrients and major elements in five catchments: (a) NH${}_{\mathrm{4}}^{+}$, (b) NO${}_{\mathrm{3}}^{-}+$ NO${}_{\mathrm{3}}^{-}$, (c) PO${}_{\mathrm{4}}^{\mathrm{3}-}$, (d) total organic N, (e) SO${}_{\mathrm{4}}^{\mathrm{2}-}$, (f) Ca²⁺, (g) K⁺, and (h) total organic C (TOC). Dashed vertical lines indicate the time of the fire.

Examining the long trends revealed that PO${}_{\mathrm{4}}^{\mathrm{3}-}$, SO${}_{\mathrm{4}}^{\mathrm{2}-}$, and K⁺ concentrations had not completely returned to pre-fire values after 3 years either in the lake (not for P) or the stream (Fig. 4). Moreover, the lake data did not show a strong response to the fire, although the stream and lake did not differ much in the pre-fire values and the whole lake catchment burned severely.

Figure 4 Long-term (13 years) changes in concentration of nutrients and major elements in Gärsjöbäcken and Märrsjön catchments: (a) NH${}_{\mathrm{4}}^{+}$, (b) NO${}_{\mathrm{2}}^{-}+$ NO${}_{\mathrm{3}}^{-}$, (c) PO${}_{\mathrm{4}}^{-\mathrm{3}}$, (d) total organic N, (e) SO${}_{\mathrm{4}}^{-\mathrm{2}}$, (f) Ca²⁺, (g) K⁺, and (h) total organic C (TOC). Dashed vertical lines indicate the time of the fire. Note that the limit of detection for K⁺ at Gärsjöbäcken was 25 µmol L⁻¹ prior to 2014.

3.3 Decay curves

Solute peaks were identified circa 1–3 months post-fire, with the two larger focus catchments (Gärsjöbäcken and Vallsjöbäcken) peaking later than the smaller catchments. Fitted solute decay curves are shown for the most intensively sampled site, Gärsjöbäcken, in Fig. 5. Summary data from the curve fitting for all four streams are shown in Table 3. For three of the four streams, the inclusion of a fast-decaying pool improved the model fits for most solutes, whereas at the strongly lake-influenced Vallsjöbäcken (flows through the largest lake), only a slow-decay pool was required to reproduce observations. Where present, the fast-decay pool contributed between 30 % and 75 % of post-fire peak concentrations, depending on site and solute, and typically had a $t\mathrm{1}/{\mathrm{2}}_{\mathrm{fast}}$ of 4–20 d. The contribution of the slow-decay pool varied very widely, from < 10 % to > 90 % of peak concentrations with a $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$ of 50–200 d.

Figure 5 Fitted decay curves for solutes exhibiting a single concentration peak after the fire in Gärsjöbäcken. See Table 3 for statistics and model fits for the other catchments.

Table 3 Decay curve modelling. Fitted initial and baseline pool concentrations, percentages, pool half-lives, and peak : baseline ratios for solutes exhibiting a single post-fire concentration peak in the four studied streams. Asterisks indicate that the 95 % credible interval of the fast-pool parameter included zero.

We observed consistent differences in the peak : baseline ratios as a function of both site and solute. In relation to site, ratios for all solutes followed the general pattern Myckelmossbäcken > Ladängsbäcken ≃ Gärsjöbäcken > Vallsjöbäcken. This sequence appears to be inversely related to the relative influence of lakes (per cent lake cover of the catchment and distance to large water body; Table 1 and Fig. 1) in the catchment upstream of the sampling point. In relation to solute, peak : baseline ratios typically followed the sequence NH${}_{\mathrm{4}}^{+}$ > SO${}_{\mathrm{4}}^{\mathrm{2}-}$ > K⁺ > TN ≥ Ca²⁺ ≃ Mg²⁺ ≃ Cl⁻. This sequence was largely replicated in the half-life data, with solutes with high peak : baseline ratios also having the shortest $t\mathrm{1}/\mathrm{2}$ values. On the other hand, we found very little evidence to suggest that $t\mathrm{1}/\mathrm{2}$ values varied consistently between the four streams.

3.4 Sensitivity analyses

Hydrological losses could have been underestimated if a flush of nutrients occurred in the first 3 weeks after the fire, prior to the start of sampling. However, the amount of precipitation was not very large in this period, so the export flux of water was low; thus solute concentration would have needed to be extremely high to generate a large solute export during this period. We consider this highly unlikely because several catchments showed solute concentration peaks a few weeks after our first sampling point, indicating that flushing (at a catchment scale) often was delayed due to buffering in the system. A sensitivity analysis for the Gärsjöbäcken catchment, assuming that the carbon and nutrient concentrations 1 week after the fire were double the values measured as the first time point, showed that the impact on the annual budget in this extreme example would nevertheless be small, resulting in an underestimation of circa 0.5 % for carbon and 1 % for nitrogen.

Treating the thin ash layer as unburned organic soil likely led to some underestimation in our carbon loss estimates due to the lower C density in ash compared to the organic soil. Using a (high) estimated ash thickness of 1 cm, a C content between 20 % and 25 %, and a wide observed ash weight (ash data from Pérez-Izquierdo et al., 2021), we calculate that treating the ash layer as unburned organic soil could have resulted in an underestimate of the average calculated carbon loss in the range of 0.01 %–1 % (2 to 45 g C m⁻²).

We did not include losses from downed wood in our C losses as this is a small component in this managed landscape. The burned area had before the fire around 4 m³ per hectare of downed wood (Jonsson et al., 2016). Assuming a stem density of 412 kg m⁻³ for Scots pine (Repola, 2006) and 50 % carbon content, the maximum loss from downed wood is on average about 80 g C m⁻² (or around 1.5 % of our calculated total C loss). This maximum value is likely an overestimation as downed wood was rarely completely consumed by the fire.

In combination, we estimate that these potential omissions in our budget calculations could have led to an underestimate of soil and forest floor total C loss of less than 3 %. Effects on budget calculations for other elements are likely smaller.

Figure 6 Conceptual C flow diagram during the first 3 years after the fire in a boreal forest catchment. Post-fire leaching is given as a range for the first 3 years, and leaching from an unburned catchment (Gärsjöbäcken, Table 2) is the mean of the 4 years before the fire occurred.

4 Discussion

4.1 Element balances

Our study shows that fire-related C and N losses resulting from a boreal wildfire were dominated by losses of the C stocks in soil O horizons, and we ascribe these losses to direct emissions during the fire (see Fig. 6 for a summary on C). Post-fire fluvial C and N losses were almost negligible compared to the deep burns in forest and peatland soils. This illustrates the importance of correctly estimating how much organic matter was consumed in the fire compared to other losses for calculating C and N budgets. The amount of C lost in the fire is around 200–1000 times higher than reported annual riverine export from boreal catchments (5–8 g m⁻² yr⁻¹; Laudon et al., 2004). We did not observe increased fluvial C losses during the first year after the fire despite the increased discharge caused by a thinner organic soil layer that decreases catchment water storage in combination with ceased plant water use. Some earlier work has suggested that fluvial dissolved C loss increases post-fire for both wildfires (Emelko et al., 2011; McEachern et al., 2000; Minshall et al., 2001) and prescribed fires (Mitchell and McDonald, 1995; Yallop et al., 2010). However, our results for TOC (considered to largely comprise DOC as discussed above) are more in line with more recent research that has found little or no effect of fire on DOC export (Betts and Jones, 2009; Burd et al., 2018; Evans et al., 2017).

Net ecosystem exchange (NEE) of CO₂ over the first 3 years post-fire indicated larger post-fire C loss than hydrologically exported C, but it still only comprised 10 % of the direct combustion emissions. It should be noted that our estimates of direct emissions may include early respiration and leaching losses, but given the severity of the fire with deep burns and large losses, we consider that combustion losses comprised the large majority of this loss. In addition, we likely underestimated direct C emissions as we did not include downed wood or biomass losses from living trees. With no vegetation, it is no surprise that the system acted as a C source immediately after fire, and the observed release of CO₂ can mainly be ascribed to heterotrophic soil respiration and to a lesser extent to dead needles and woody biomass. Compared to undisturbed systems, heterotrophic respiration actually seems to decrease after fire (reviewed in Amiro et al., 2003) partly due to the formation of inert carbon, i.e. pyrogenic carbon that may stabilize the remaining organic carbon (Jones et al., 2019). About 3 years post-fire, summer NEE showed for the first time net C uptake. It is likely that the overall pattern was similar across the whole burn because we observed a rapid increase in LAI in all catchments. Flux data from boreal North America have also shown summer net C uptake 2 years post-fire, but it may take 10 years until the system is a sink on an annual basis (Amiro et al., 2003, 2010; Goulden et al., 2011; Kashian et al., 2013). Specifically, an eddy covariance study in boreal Canada estimated the net ecosystem production 1 and 2 years post-fire and reported C losses of 192 and 93 g C m⁻² yr⁻¹, respectively (Goulden et al., 2011). These values are similar to our two sites (155 to 165 g C m⁻² yr⁻¹ over 2 years), but further research is needed to establish if such values are typical of boreal uplands post-fire. However, despite our effort to track carbon flows in the system, we still had to model flux values for the first fall–winter period, and combustion losses were inferred by using unburned reference plots. Hence, our estimates are associated with uncertainty that needs to be considered when upscaling these results.

In contrast to C, we observed a dramatic increase in hydrological N loss that was largely driven by higher concentrations in the streams. The amount of dissolved N lost over the first years (almost 1 g N m⁻²) may be small compared to the direct combustion losses (< 1 %), but this is available N, whereas much of the N lost in the fire is N which would have been derived from forms of stable organic matter that were not readily available for the plants (Smith et al., 2011; Tamm, 1991). Our annual estimates of fluvial N losses are similar to those reported for a mixed coniferous forest (Nevada, US; Johnson et al., 2007) and for peaty heathland (Northern Ireland, UK; Evans et al., 2017), but there are losses 100 times greater than had been reported for a Mediterranean shrubland (Dannenmann et al., 2018). Our estimates of direct N losses are at the higher end of reported values for temperate and boreal coniferous forests, 30 to 90 g m⁻² (Brais et al., 2000; Grier, 1975; Johnson et al., 2007), but in contrast to previous studies, we included N losses from drained peatlands that probably resulted in higher total losses. In addition, both other estimates and our own do not include post-fire gaseous N emissions that during the first post-fire year have been shown to comprise 10 %–15 % of the direct fire combustion losses in shrubland systems (Dannenmann et al., 2018).

Despite these large N losses, there is little evidence that either direct or fluvial N losses are relevant for post-fire productivity at a catchment scale. Recently, a study by Turner et al. (2019) showed a remarkably rapid post-fire (4 years) build-up of soil N and little evidence that the N loss had a long-term impact on productivity. It is hypothesized that post-fire plant communities, if quickly established, can retain N before it is lost hydrologically (Smithwick et al., 2009). At our study site, vegetation established after 2 years, but most of the soluble N had already been lost by that time. Hence, plants must utilize newly mineralized N or acquire their N through microbes (e.g. via N-fixation). In fact, it is unknown how plants can acquire large amounts of N post-fire and how the N pool builds up quicker than estimated N-fixation rates (Turner et al., 2019). Our N losses (ca. 100 g m⁻²), for example, correspond to more than 150 years of N input from fixation and deposition (based on 0.6 g m⁻² yr⁻¹ N input) (Brais et al., 2000; Zackrisson et al., 2004). Clearly, fire is a key driver of the global N cycle. Future studies should focus on elucidating the mechanisms behind post-fire N build-up in the boreal biome to better capture this dynamic in ecosystem models.

Few studies have quantified other fire-related nutrient losses such as S, P, K, and Mg. The integrated hydrological mass export during the first year after the fire corresponds to around 5 years (P, K, and Mg) and 26 years (S) of pre-fire element export (Table 2). Hence, on a longer timescale, these losses seem unlikely to affect the productivity of the system, although they could influence short-term availability for uptake by the biota, as well as soil acidity, in these relatively base-poor ecosystems. Instead, our study indicates that soil and biomass retention capacity for base cations was fast and efficient in this fire-impacted boreal ecosystem. Shorter fire intervals might therefore have a limited impact on base cation budgets, although it is clear that they will fundamentally alter C and N budgets due to loss of slow-forming organic soil.

4.2 Water quality and decay curves

The fire had generally a strong short-term impact on the water quality with large short-term variations of both base cations and acid anions over time. Hence, our study highlights the importance of frequent sampling soon after the fire to accurately capture the post-fire dynamics in water chemistry. At these peatland-rich sites, pH remained fairly stable despite the great fluctuations in mineral anions (SO${}_{\mathrm{4}}^{-\mathrm{2}}$, NO${}_{\mathrm{3}}^{-}$; Fig. 3). A significant pH drop only occurred at the peatland-dominated site of Myckelmossbäcken where TOC was initially suppressed just after the fire (Figs. 3, S3). At all other sites organic anion concentrations were above 100 µeq L⁻¹ (Fig. S4) which buffered pH against any potential charge imbalance of sulfate and base cations. TOC is mainly released from riparian peatlands in boreal catchments (Ledesma et al., 2015), and it is possible that an intact (less burned) riparian zone through its TOC release can buffer and thus prevent a large pH drop from occurring after fire.

Sustained elevated levels of reactive phosphorus have been reported for other boreal wildfires, and our relative increase are similar to studies examining phosphorus concentration up to 5 years post-fire (Hauer and Spencer, 1998; Silins et al., 2014). The P and N enrichment likely caused higher algal productivity in streams, which can generate effects at higher trophic levels (Silins et al., 2014), but this was not monitored in our study.

The analysis of decay curves suggests that there are two distinct sources of solute flushing to the stream. A single-exponential model was unable to reproduce both the rapid initial decline and the longer-term decrease, whereas a two-pool model generally gave a good fit. Moreover, a two-pool model is mechanistically interpretable. The first “fast-decay” pool is associated with the immediate post-fire period, typically made a significant contribution to peak solute concentrations, and was observed for most solutes in three of the four streams, with the exception of the strongly lake-influenced Vallsjöbäcken. The rapidity with which this peak dissipated, with half-lives between 4 and 25 d, suggests that it reflects the instantaneous mobilization of solutes due to pyrolysis of biomass and soil organic matter, followed by hydrologically controlled flushing into the drainage network. The fine ash that formed is most probably very soluble and may be leached out fast with rainwater (Grier, 1975). The second, “slow-decay” pool contributed variably to post-fire peak concentrations but affected water chemistry for a period of years, with half-lives typically in the order of 75–175 d. The consistent differences in $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$ between solutes, coupled with the absence of clear variability in $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$ between streams, leads us to conclude that this pool is largely determined by biogeochemical processes occurring after the fire. This fits with the observed heterotrophic respiration in our NEE data and suggests gradual leaching of solutes from ash and the breakdown and dissolution of dead organic matter. The relative contribution of the two pools of element leaching is likely determined by burn severity, in which a more severe burn would increase the size of the fast pool by consuming more of the organic matter, leaving the inorganics (K, Ca, NH${}_{\mathrm{4}}^{+}$, etc.) available for rapid leaching. This would also suggest that more severe fires result in a smaller “slow” pool because there is less organic matter left to decompose.

Differences in peak : baseline ratios and $t\mathrm{1}/{\mathrm{2}}_{\mathrm{slow}}$ between solutes appear to reflect their source within the ecosystem; N and K are largely present in non-woody biomass including microbes, leaves, and fine roots and are therefore likely to be released relatively quickly. In particular, NH${}_{\mathrm{4}}^{+}$ is the initial product of organic matter mineralization, and the very large and fast-declining peaks observed in this solute (e.g. compared to either NO${}_{\mathrm{3}}^{-}$ or TN) suggest that the supply of NH${}_{\mathrm{4}}^{+}$ due to organic matter pyrolysis and mineralization immediately following the fire overwhelmed abiotic and biotic retention mechanisms, as well as terrestrial and aquatic nitrification capacity. This short-lived NH${}_{\mathrm{4}}^{+}$ pulse, together with more sustained leaching of NO${}_{\mathrm{3}}^{-}$ in the years after the fire, is consistent with previous studies of wildfire impacts (e.g. Wan et al., 2001) and with other studies of N cycle responses to major ecosystem disturbances, such as bark beetle attacks (Kopáček et al., 2018). The mechanisms behind such similar responses to different disturbances are likely less plant uptake and increased N mineralization. In contrast to N solutes, the divalent base cations are more structurally bound within biomass pools, strongly retained on soil cation exchange sites, and therefore released more gradually via organic matter mineralization, especially in the presence of pyrogenic organic matter. The slow release of Cl also suggests release from decaying organic matter, consistent with previous studies suggesting that large amounts of Cl is biotically cycled within northern forest ecosystems (Bastviken et al., 2006). The source of SO${}_{\mathrm{4}}^{-\mathrm{2}}$ leaching may be somewhat different because the largest pools of S in our study catchments are believed to be sulfides and organic S compounds held under anaerobic conditions in wetlands (Schiff et al., 2005). Thus the largest peaks in SO${}_{\mathrm{4}}^{-\mathrm{2}}$ were recorded in the peat-influenced Myckelmossbäcken and are associated with the combustion of a considerable depth of peat.

In two of the investigated stream catchments most of the forest stands were salvage logged during the first year after the fire. Interestingly, we did not observe any clear or consistent differences in water quality between salvage-logged and non-salvage-logged catchments over the study period. A study by Silins et al. (2014), possibly the only study that has made this comparison for boreal catchments, found larger increases in stream P concentration in salvage-logged catchments. However, this was in an area with extreme topography (Rocky Mountains) where mechanical damage led to increased erosion. In our lower-relief study area, evidence of large-scale soil disturbance during salvage logging was not observed. We did not investigate post-fire CO₂ fluxes in logged areas, but previous studies have not found clear evidence of increased soil C losses compared to unlogged areas (Kishchuk et al., 2016; Parro et al., 2019). To better investigate if post-fire salvage logging has an ecologically important effect on water quality in boreal Europe more catchments, and longer time series are needed.

5 Conclusions

Our study provides a unique integrated quantification of the impact of wildfires on boreal forest biogeochemistry (e.g. Fig. 6). Overall, hydrological export of nutrients was fairly short-lived (1–2 years) and was caused mainly by higher ion concentrations and not by increased discharge. For some major elements with gaseous loss pathways, notably C and N, fluvial losses were small compared to the direct emission; in fact, no increase in aqueous C export was observed, and fluvial losses of N and C can be considered minor compared to combustion losses for boreal catchment budgets during a fire. For elements that showed elevated exports (N, P, S, Mg, K), the first year post-fire was equivalent to circa 5 years (26 for S) of exports in unburned systems. Base cation fluxes 3 years post-fire were similar to pre-fire conditions except for K that remained elevated much longer, suggesting slower release and weaker retention of this element. Our decay curves and comparable pre- and post-fire fluxes indicate that the boreal forest ecosystem has re-established a similar steady-state of deposition, weathering, and export. It will be interesting to revisit these catchments in a few years to study whether the element uptake of the growing trees will lead to lower stream water export. Vegetation regrowth was rapid and likely contributed to decreased leaching of nutrients while initiating C sequestration of the system. After 3 years post-fire, there was a clear net ecosystem C uptake during the summer, suggesting that fire-induced C losses had largely concluded and that the ecosystem will likely become a net CO₂ sink in future years as the forest regrows. However, given the magnitude of C loss from the combustion of the organic soil, it will likely take decades or even centuries for overall ecosystem C stocks to recover. If fire frequency increases across boreal forest ecosystems, these forests can become net long-term sources of CO₂ to the atmosphere, reversing their current function as carbon sinks.