Introduction

The role of forests in regulating greenhouse gas (GHG)1 budget, in particular for nitrous oxide (N2O), is still largely unknown2. The increase in atmospheric N2O levels from a pre-industrial concentration of 270 ppbv to 328 ppbv in 20163 is of concern because: (i) N2O is responsible for approximately 6% of global radiative forcing from anthropogenic GHGs, (ii) N2O is the ruling stratospheric ozone-layer depleting agent in the 21st century4, and (iii) N2O the third most important GHG with a global warming potential 296 times that of CO23 (100-year lifetime adjustment, with feedbacks). A variety of nitrogen cycle processes can produce N2O5 but in riparian soils, denitrification is the most important source of N2O6.

Hot spots and hot moments (extreme emission events) largely determine spatio-temporal variation of N2O emissions from soils7, and soil volumetric water content (SWC) is a leading factor controlling hot moments. An SWC value of 0.5–0.8 m3 m−3 has been shown to be optimal for soil N2O emissions in forests on both mineral soils8 and organic soils9. Therefore, depending on the initial moisture, both flooding and drought can induce hot moments in forests10. Drought decreases soil N2O emissions or even turns the soil into a sink of N2O11. Another mechanism that creates hot moments of soil N2O emission is freeze-thaw cycles12,13,14,15. Although there are several hypotheses explaining the impact of freeze-thaw cycles on soil N2O emissions14,16,17, a generally accepted theory of the impact of freeze-thaw on N2O fluxes is still missing.

Besides soils, plants including trees are known to exchange N2O with the atmosphere18,19,20,21. N2O emitted from tree stems and canopy may originate from the soil18,19, from microbial or fungal activity within the tree stem, or from plant biophysical processes21. The relative contributions of these sources to N2O exchange between stems and the atmosphere is unclear in any tree species, including Gray alder (Alnus incana), a common species in riparian areas that are widespread across the Northern Hemisphere22 (Supplementary Fig. 1).

Above tree canopies, long-term ecosystem-level N2O flux measurements using the eddy covariance (EC) technology are still scarce23,24,25. No previous complex investigations on forest ecosystems’ N2O budgets (soil and tree stem flux studies with EC measurements above the canopy) exist.

We analyze 2.5-years seasonal and interannual dynamics of N2O fluxes in a deciduous riparian forest using high temporal frequency measurements made at the scale of the whole ecosystem and from two ecosystem components, soils and tree stems. We relate the fluxes to key environmental factors to test the following hypotheses: (1) soil water content (SWC)- and temperature-related hot moments determine long-term patterns of N2O fluxes, and (2) EC-measured N2O fluxes are coherent with the sum of soil and stem fluxes.

Results

Soil N2O fluxes

Across the study period, the soil N2O fluxes were highly dynamic, with the annual budget largely driven by the short-lived hot moments (Fig. 1 and Supplementary Fig. 2). N2O fluxes varied from −0.040 to 1.50 mg N m−2 h−1 while daily averages reached 0.225 mg N m−2 h−1 (Fig. 2). The heatmap in Supplementary Fig. 3 presents a spatial and temporal variation of these values, showing that across the whole study period, no remarkable difference between the values was measured in individual chambers. Negative N2O fluxes occurred primarily in winter months when the near-ground air temperature was consistently negative, accounting for 43% of monthly values in February 2018 and January 2019 (Supplementary Fig. 4). We identified four hot moments in soil N2O flux time series defined by remarkably higher flux values (Fig. 1). The highest flux values (>0.10 mg N2O–N m−2 h−1) were observed during the Drought Onset period starting in May 2018, accounting for a maximum of 38% of the measurements (Supplementary Fig. 3). During this short period, the daily average flux in all individual chambers was >0.02 mg N m−2 h−1 (Supplementary Fig. 3). The freeze-thaw period in 2019 contributed 14% of the whole soil emission (Supplementary Table 1). We attribute the absence of a freeze-thaw hot moment in 2018 to consistently negative near-ground air temperature and stable (2–12 cm deep) snow cover in January and February 2018 (Fig. 1) that would impede N2O diffusion rates.

Fig. 1: Dynamics of ecosystem-level N2O fluxes in the Agali gray alder forest during the study period Sep 2017–Dec 2019.
figure 1

Lines are 5-day median values and shaded areas are 25th and 75th percentiles. Vertical lines show the beginning and end of hot moments: 1 = Wet (2017-09-01 … 2017-12-01), 2 = Dry (2018-05-01 … 2018-08-05) with 2a = Drought Onset (2018-05-02 … 2018-05-21), 3 = freeze-thaw (2019-02-01 … 2019-03-15), and 4 = Dry-Minor (2019-06-01 … 2019-06-30). The air temperature was measured at 50 cm from the soil surface representing the near-ground air temperature. Note that the origins are on opposite ends of the left- and right-side axes.

Full size image
Fig. 2: Relationships between the soil temperature, soil water content (SWC), and flux of N2O from the soil over the study period.
figure 2

a Contour plot showing relationships between soil temperature, SWC, and N2O emission (n = 755). b Regression curve of SWC vs. N2O fluxes. Curve fitted regression of SWC and N2O flux (R2 = 0.07, p < 0.01, n = 757). N2O = (15725.05 × SWC7.73) × exp(−15.38 × SWC). c Regression curve of soil temperature vs N2O. Local polynomial regression fitting of soil temperature and N2O flux (R2 = 0.13, p < 0.01, n = 756). The dashed red lines represent 95% confidence intervals for the regression line.

Full size image

During the whole study period (September 2017–December 2019) the cumulative N2O soil flux was 458.8 ± 7.7 mg N m−2 (mean ± standard error, SE) and the hot moments contributed 56% of the whole flux (Supplementary Table 1). During the calendar years 2018 and 2019, 196.3 ± 7.1 and 221.0 ± 12.4 mg N m−2 year−1 were emitted as N2O (Supplementary Table 1). Considering the two full years of the study (September 2017– September 2018 and September 2018–September 2019), the corresponding cumulative flux values were 215.5 ± 7.7 and 221.4 ± 12.2 mg N m−2 year−1 (Fig. 3; Supplementary Table 2)

Fig. 3: Cumulative fluxes of N2O from the soil, stems, and the ecosystem (eddy covariance above the canopies).
figure 3

Data from two full years (September 2017–September 2019) and one half-year (September–December 2019). Due to significantly lower values, the stem fluxes are also plotted in the inset. Notice that the stem fluxes have been measured from Sept 2017 to Dec 2018.

Full size image

Except for the dry hot moments, we found no remarkable diurnal pattern of soil N2O fluxes. Average day-time values during the dry hot moment were up to 100 times higher than those in the night-time showing that the variability was higher than in the other months (Supplementary Fig. 5).

Stem N2O fluxes

Stem fluxes of N2O measured during 52 campaigns averaged over all measured heights and expressed per m2 of soil surface varied between −0.00028 and 0.0228 mg N m−2 h−1. The highest emissions were measured on the lowest position of tree stems whereas at higher positions (170 cm aboveground) slight consumption was observed. Stem flux during the measurement period from September 2017 to December 2018 was 0.00022 ± 0.00007 mg N2O–N m−2 h−1 (mean ± SE). The wet and dry hot moments contributed respectively 41% and 12% of total cumulative emissions during the whole stem measurement period (3.53 mg N m−2; Supplementary Table 1; Fig. 3).

Ecosystem level N2O fluxes

The daily sums of ecosystem N2O fluxes varied relatively little, ranging from −0.60 to 1.16 mg N m−2 d−1. Of the whole measurement period, 76% passed the quality check and we gap-filled the remaining 24% (Fig. 4). Hot moments of soil N2O fluxes were not reflected at the ecosystem scale (Figs. 1 and 3), but showed a seasonal pattern of high positive fluxes in spring (March–April) and autumn months (October–November) and small fluxes near zero over summer months (Figs. 4 and S6). Likewise, higher values were observed during the “Wet” hot moment with 22% of the total cumulative EC flux of the whole study period (78.2 mg N m−2; Figs. 3 and 4, Supplementary Table 1). We observed a distinctive diurnal pattern with small negative fluxes during the morning hours (8 a.m. to 12 a.m.) in the summer months of both years. We observed no diurnal pattern in either the simultaneous soil fluxes or in EC fluxes of autumn, winter, or spring months. The hourly average EC flux values ranged from −0.05 to 0.03 mg N2O–N m−2 h−1 (Supplementary Fig. 4).

Fig. 4: Seasonal cycle of ecosystem N2O flux measured with quantum cascade laser absorption spectrometer in eddy tower.
figure 4

The markers denote daily total values, the line is a 7-day running mean. The periods marked with red color represent time intervals with gap-filled data (MDS-method) exceeding 50%. Hot moments: 1—wet, 2—drought (without showing drought onset), 3—freeze-thaw, and 4—dry-minor.

Full size image

Relationships between N2O flux and environmental parameters

The main factors related to N2O soil fluxes in this riparian forest ecosystem were SWC and soil temperature (Fig. 2a). Based on the full data set measured during the study period, there was an optimal SWC value of about 0.5 m3 m−3 at which the soil flux was the highest (Fig. 2b). The relationship between the soil surface temperature and soil N2O flux was more complex showing two peaks, one at 0–4 °C and a second one at 13–14 °C (Fig. 2c). The first peak corresponds to the freeze-thaw hot moment and the second one represents the dry and dry-minor hot moments.

During the Dry hot moment, the correlation between the SWC and soil N2O emission was very strong, showing a clear peak at SWC values between 0.35 and 0.5 m3 m−3 (Fig. 5).

Fig. 5: Dynamics of soil water content (SWC) vs soil N2O flux (hourly average values) during the dry hot moment (2018-05-01 … 2018-08-05).
figure 5

The curve is calculated after the Bragg equation with four parameters: Y = c + (d − c) × exp(−b)× (X − e)2, where b = 92.77, c = 0.0123, d = 0.1786, e = 0.4314 and X is SWC and Y is N2O flux (R2 = 0.74, p < 0.001, n = 1065).

Full size image

We did not find any significant relationships between ecosystem N2O fluxes and environmental parameters (air temperature, soil temperature, SWC, wind speed) or gross primary production (GPP) on the half-hourly scale. However, weekly average fluxes followed the same general pattern as that of SWC with the modifying influence of changes in air temperature (Supplementary Fig. 6). The stem N2O fluxes were also related to SWC with a positive but insignificant correlation (p = 0.09) while no clear relationship was found during the hot moments.

Comparing the monthly average soil and EC fluxes of N2O to weather patterns, we found that during the spring (drought onset and end of freeze-thaw periods) the largest differences between the soil and EC fluxes coincided with the largest percentage of clear days (cloudiness ratio <0.4; Fig. 6). Differences between the soil and EC fluxes were lowest during the autumn and winter with the highest number of days of potential fog formation (dew point depression (DPD) < 4 °C) (Fig. 6).

Fig. 6: N2O fluxes (soil and eddy covariance) related to the chance of sunshine and fog formation.
figure 6

Relationship between the monthly sum of soil and EC fluxes of N2O (mg N m−2 month−1) with the percentage of days with a high chance of sunshine (cloudiness ratio <0.4) and fog formation (dew point depression <4 °C). During the drought onset (May 2018) and dry-minor (June–July 2019) periods photochemical reactions might play an important role in decreasing ecosystem (EC) flux of N2O, whereas during the Wet (September–November 2017) and freeze-thaw (February 2019) periods absorption on wet surfaces and water droplets could be the reason for the decrease in EC N2O flux. The difference between the soil and EC N2O fluxes is negatively correlated with the percentage of days in a month with dew point depression <4 °C, whereas a positive trend is observed both in the difference between soil and EC flux and the share of sunny days (cloudiness ratio < 0.4) in a month.

Full size image

Discussion

We hypothesized that SWC was the main factor associated with hot moments of N2O emission from a riparian forest ecosystem. Our data generally supported this hypothesis for soil emissions, but we suggest that drought and rewetting triggered these events by different mechanisms as previously argued26. We identified four periods during which N2O hot moments arose due to distinct mechanisms (Figs. 1 and 2). Across all cases, N2O fluxes peaked at SWC value around 0.5 m3 m−3.

Soil N2O emissions during the dry and freeze-thaw hot moments were 8–10 times higher than the period averages (Supplementary Fig. 2). N2O emissions from intact soil mesocosms from a temperate forest showed an increase of up to four orders of magnitude following flooding pulses27, consistent with the effects of fluctuating groundwater table depth on N2O emissions under field conditions28.

The highest soil N2O fluxes in the present study were observed at SWC values around 0.5 m3 m−3 (Fig. 2b), a pattern that was most dramatic during drought onset (Fig. 5). Although several studies have shown a similar relationship8,9, the direction of change in soil N2O emissions depends on whether the soil is initially dry (SWC < 0.45 m3 m−3) or wet (SWC > 0.45 m3 m−3). Drought has been shown to decrease N2O emissions in dry mineral soils15, but increase N2O emissions in wet soils (SWC > 0.6 m3 m−3)9,29. We observed both trends: a substantial drought-driven increase of N2O emissions under wet conditions (SWC > 0.7 m3 m−3; Drought Onset episode) and a drought-driven decrease of N2O emissions in dry conditions (SWC < 0.45 m3 m−3; end of dry period) with short-term emission peaks caused by precipitation-driven rewetting (dry-minor period; Supplementary Fig. 2). In all cases, the highest emissions occurred when the SWC transitioned past an optimum value of about 0.5 m3 m−3 (Supplementary Fig. 2).

During the freeze-thaw period, a different complex set of factors caused a substantial increase in soil N2O emissions. Several hypotheses have been posed to explain this phenomenon, the most common of which are: (i) disruption of soil aggregates exposing physically protected organic matter to rapid mineralization by microorganisms14; (ii) the death of microorganisms, fine roots, and mycorrhiza during the freeze16, providing rapidly decomposable organic matter for the thaw; and (iii) the death of fine roots decreasing competition between roots and microorganisms for nitrate, the main substrate for N2O production16. However, the underlying mechanism involved in the pulse emissions of N2O remains uncertain and further exploration is required30. In addition, rising SWC played a role in N2O fluxes at the end of the freeze-thaw (Supplementary Fig. 3), similarly to the observations of Öquist et al.31.

The different hot moment mechanisms were driven by different relationships with the environmental factors. During the drought onset, there was a strong negative correlation (RSpearman = −0.93) between the speed of SWC decrease and N2O flux (Supplementary Fig. 7a). Only a few authors have considered drying speed as a factor in soil nitrogen cycle processes32,33,34 and we report the phenomenon in situ. Accordingly, the speed of SWC decrease could be included in N2O flux models. This resonates with the review by Barrat et al.26 showing that the bigger the increase in SWC the larger the N2O pulse. In our freeze-thaw period, a positive correlation (RSpearman = 0.64) between the near-floor air temperature and N2O flux was observed (Supplementary Fig. 7b) suggesting that temperature may play a dominant role in the initiation of the hot moments. Likely, after the freeze, which physically breaks soil organic matter, and releases nutrients and dead microbial biomass, an increase in temperature stimulate microbial activity. This is fed by the newly released carbon and nitrogen to form N2O peaks. Indeed, we did not observe meaningful diurnal patterns in soil N2O fluxes during most of the study period, except for the dry and freeze-thaw hot moments (Supplementary Fig. 4). Until now only a few studies on terrestrial ecosystems have reported differences in soil N2O emissions between day and night35.

The second hypothesis that EC-measured N2O fluxes are coherent with the sum of soil and stem fluxes was not supported. The ecosystem N2O flux measured by the EC technique above the forest canopy was relatively low and did not follow the variability in soil N2O fluxes (Figs. 1, 2, and 4). During the wet period, one significant N2O peak was observed, while emissions remained low during the rest of the period. Other hot moments of soil N2O flux did not increase ecosystem N2O flux (Fig. 4, Supplementary Figs. 2 and 8). The cumulative N2O emission from the ecosystem (78.2 mg N m−2) was 5.9 times smaller than the soil’s emission (458.8 mg N m−2; Supplementary Table 1; Fig. 3). During the 1.5 years, the cumulative flux from alder stems was 3.53 mg N m–2, which constituted only 0.8% of cumulative soil fluxes (Supplementary Table 1).

Because the measurement frequency was sufficient to capture all the fluxes above the canopy it is difficult to explain such a large difference between the soil and ecosystem emissions. We offer the following explanations of this difference:

  1. (1)

    Horizontal advective fluxes and decoupling of gas fluxes between the layers are likely important reasons. While our study site is located on flat terrain, below-canopy horizontal advection could still play a role, especially in a closed canopy forest ecosystem. Dense canopies can reduce vertical mixing and decouple the layers below and above the canopy36 in a way that peaks of soil N2O emissions leak horizontally undetected by the EC system. For instance, during the Drought Onset episode in 2019, N2O mixed ratios increased during wind-still nights (Supplementary Fig. 8) undetected by the EC measurements. Probably, a combination of poor mixing of soil-borne N2O during low friction velocity (u*) and wind speed and decoupling of soil-borne and above canopy N2O flux measurements due to the dense canopy layer (Supplementary Fig. 8). Nevertheless, the absence of concentration profile measurements and our simplified approach to the storage calculation may have contributed to the discrepancy between the soil chamber and EC measurements.

  2. (2)

    UV-induced photodissociation may play a role in decreasing N2O flux above the canopy. Typically for the drought onset episode, EC fluxes declined during the clear sunny days (Supplementary Fig. 4b). Thus, reduction during daylight could be explained by the UV-induced photodissociation of N2O as found in other ecosystems37. Indeed, UV-induced N2O photolysis is included as an important process in the assessment of the atmospheric lifetime of N2O38. In our study, the monthly sum of sunshine hours was highest during the drought onset and dry-minor episodes when the differences between the soil and EC fluxes of N2O were large (Fig. 6; Supplementary Fig. 9). During the dry period from May to August 2018, midday UV-B radiation of 295–315 nm wavelength was consistently 1–1.5 W m−2 (unpublished data from the nearby SMEAR Estonia station), up to 6 times higher than the cloudless noon values measured in Scandinavia by Fraser et al.39 and in 2016 in our study area40. Thus, during the dry and dry-minor episodes, photochemical reactions may have decreased the ecosystem (EC) flux of N2O measured above the canopy. On the other hand, the canopy shelter could mitigate the effect of UV radiation and increase the N2O concentration within the canopy. This again proves the need for profile measurements in such studies.

  3. (3)

    Potential N2O dissolution in atmospheric water. In autumn, winter, and spring the consumption of N2O could be related to the high solubility of N2O gas in water (1.0 ml gas per ml water at 5 °C)41. For example, Warneke et al.42 reported that absorption of N2O in woodland soil water contributed up to half of the total N2O consumption in the soil. In the present case, part of the N2O produced during the freeze-thaw period may have been absorbed in fog. Eugster et al.24 in their early EC study in a mixed forest suggested that wetting of the canopy in fog can have a strong influence on N2O fluxes; however, no clear evidence of absorption in fog has been reported. Likewise, Min et al.43 found that NOx fluxes from the forest canopy were smaller than measured soil NOx emissions and referred to the phenomenon as a “canopy reduction factor” which they applied to soil NOx emissions in large-scale models. The interpretation of these differences was a chemical conversion of NOx to other nitrogen oxides within the forest canopy. Fulgham et al.44 report that wet surfaces of leaves, needles, and branches in a mixed forest control the vertical exchange of volatile organic acids. The exchange velocity of these gases was well correlated with DPD. This is compatible with our findings on the difference between the soil and EC fluxes during the freeze-thaw period with high DPD (Fig. 6; Supplementary Figs. 9 and 10). We speculate that some of the soil-emitted N2O is absorbed in the water droplets and films on soil, stem, and leaf surfaces, hence removing N2O from the air. From the canopy, dissolved N2O can be transported to soils via stemflow and throughfall45. We also observed either correspondence or a couple-week difference between the dynamics of N2O EC flux and the peaks of its possible regulator—the maximum potential concentration of dissolved N2O in m3 of air calculated based on LiCor-measured water vapor at EC tower (Supplementary Fig. 11).

It is difficult to evaluate the relative importance of these non-exclusive causes of the discrepancy between the soil, stem, and ecosystem fluxes of N2O with the present data set. In addition to the soil and tree stem measurements made in this riparian forest, a better understanding of forest N2O budgets also require the canopy-, and leaf-level measurements as well as assessment of advective flux and storage by profile measurements24,46.

Comparison of our results with other studies is hampered by the fact that almost all previously reported annual rates of N2O emissions from forests are based on soil surface measurements only. Our alder soil surface emitted an average of 2.18 kg N ha−1 year−1 (average for the period Sept 2017–Sept 2019; Supplementary Table 2) similar to the IPCC emission factors for the boreal drained nutrient-rich forest on organic soil47 (3.2 kg N ha−1 year−1) and temperate drained forests on organic soil47 (the 2.8 kg N ha−1 year−1), but several times less than from peatlands drained for agriculture48 or other drained N-rich peatlands9. However, based on our whole-ecosystem EC data we found that above-canopy N2O emissions are an order of magnitude less (0.35 kg N ha−1 year−1; Supplementary Table 2) than the soil estimates. This implies an 84% overestimate of soil-based emission estimates for temperate drained forests on organic soils.

Gray alder (Alnus incana) forests are widely distributed in Europe and North America22 (Supplementary Fig. 1), often dominating riparian forest communities49. Upscaling our N2O flux values to the whole Alnus incana subsp incana distribution area in Europe (15,000 km2)22, we estimate total annual emissions of 3270 (soil-based) or 390 tons (forest ecosystem-based) of N a year. Thus, in addition to several ecosystem services which riparian alder forests can provide, they are low emitters of N2O and therefore attractive for riparian zone management when considering climate impacts.

The outcomes of our 2.5-years study show that this alder forest is a weak net source of N2O. Hot moments related to SWC and temperature accounted for 56% of soil emissions throughout the whole study period. The soil N2O flux peaked at about 0.5 m3 m−3 of SWC. The ecosystem (EC) flux was not coherent with the sum of soil and stem flux. In comparison to the high soil N2O emission, the ecosystem level emission was about 5.9 times lower. Photochemical reactions and dissolution in atmospheric water are potential mechanisms by which soil-emitted N2O is consumed in the forest canopies. Our results suggest that the reported N2O emissions from temperate drained forests on organic soils have been greatly overestimated.

In the next decades, we anticipate a global increase in the frequency of disturbances causing hot moments of GHG emissions in terrestrial ecosystems. Our study reveals the importance of high-frequency ecosystem-scale field measurements across the year in order to accurately quantify forest GHG emissions. Full understanding of nitrogen budgets of riparian forests cannot rely only on soil level measurements but must be combined with tree-stem, canopy, and ecosystem-level EC fluxes.

Methods

Study site and set-up

The studied hemiboreal riparian forest is a 40-year old Filipendula type gray alder (Alnus incana (L.) Moench) forest stand grown on a former agricultural land. It is situated in the Agali Village (58°17′N; 27°17′E) in eastern Estonia within the Lake Peipsi Lowland50 (Supplementary Figs. 12 and 13).

The area is characterized by a flat relief with an average elevation of 32 m a.s.l., formed from the bottom of former periglacial lake systems, it is slightly inclined (1%) towards a tributary of the Kalli River. The soil is Gleyic Luvisol. The thickness of the humus layer was 15–20 cm. The content of total carbon (TC), total nitrogen (TN), nitrate (NO3−–N), ammonia NH4+–N, Ca, and Mg per dry matter in 10 cm topsoil was 3.8 and 0.33%, and 2.42, 2.89, 1487 and 283 mg kg−1, respectively, which was correspondingly 6.3, 8.3, 4.4, 3.6, 2.3, and 2.0 times more than those in 20 cm deep zone (Supplementary Table 3).

The long-term average annual precipitation of the region is 650 mm, and the average temperature is 17.0 °C in July and –6.7 °C in January. The duration of the growing season is typically 175–180 days from mid-April to October51.

The mean height of the forest stand is 17.5 m, stand density 1520 trees per ha, the mean stem diameter at breast height 15.6 cm, basal area 30.5 m2 ha−1, the growing stock 245 m3 ha−1, and the current annual increment of stems 12.0 m3 ha−1 year−1 (based on Uri et al.52 and Becker et al.53). In the forest floor, the following herbs dominate: Filipendula ulmaria (L.) Maxim., Aegopodium podagraria L., Cirsium oleraceum (L.) Scop., Geum rivale L., Crepis paludosa (L.) Moench, shrubs (Rubus idaeus L., Frangula alnus L., Daphne mezereum L.), and young trees (A. incana, Prunus padus (L.)) dominate. In the moss-layer Climacium dendroides (Hedw.) F. Weber & D. Mohr, Plagiomnium spp and Rhytidiadelphus triquetrus (Hedw.) Warnst are overwhelming.

Environmental characteristics of hot moments

Based on high emissions of N2O, dynamics of SWC, and near-ground air temperature, we identified four hot moments and related them to soil and environmental variables (see numbers in Fig. 1): wet (1), dry (2) with drought onset (2a), freeze-thaw (3), and dry-minor (4). The main criterion for the hot moments was a rapid increase in N2O emissions of any source.

Anomalies from the mean of each hot moment period illustrate the pattern of fluxes during the hot moments (Supplementary Fig. 2). At the end of the freeze-thaw period, the rising SWC is driven by snowmelt became a leading determinant (Supplementary Fig. 2). During the wet period, the rise in soil emissions was accompanied by a remarkable increase in the EC-based ecosystem fluxes. However, all the other hot moments were isolated to soil surfaces.

Soil flux measurements

Soil fluxes were measured using 12 automatic dynamic chambers located at 1–2 m distance from each studied tree and installed in June 2017 (Supplementary Fig. 11, see also54). The chambers were made from polymethyl methacrylate (Plexiglas) covered with non-transparent plastic film. Each soil chamber (volume of 0.032 m³) covered a 0.16 m² soil surface. To avoid stratification of gas inside the chamber, air with a constant flow rate of 1.8 L min−1 was circulated within a closed loop between the chamber and gas analyzer unit during the measurements by a diaphragm pump. The air sample was taken from the top of the chamber headspace and pumped back by distributing it to each side of the chamber. For the measurements, the soil chambers were closed automatically for 9 min each. The flushing time of the whole system with ambient air between measurement periods was 1 min. Thus, there were ~12 measurements per chamber per day, making a total of 144 flux measurements per day. A Picarro G2508 (Picarro Inc., Santa Clara, CA, USA) gas analyzer using cavity ring-down spectroscopy (CRDS) technology was used to monitor N2O gas concentrations in the frequency of ~1.17 measurements per second. The chambers were connected to the gas analyzer using a multiplexer allowing a sequent practically continuous measurement.

To account for initial stabilization after chamber closing and flushing time, we used 5 min out of the total 9 min closing time (~350 concentration measurements) to estimate slope change of N2O concentration, which was the basis for soil flux calculations.

After the quality check, 105,830 flux values (98.7% of total possible) of soil N2O fluxes could be used during the whole study period.

Stem flux measurements

The tree stem fluxes were measured manually with frequency 1–2 times per week from September 2017 until December 2018. Twelve representative mature gray alder trees were selected for stem flux measurements and equipped with static closed tree stem chamber systems for stem flux measurements20. Soil fluxes were investigated close to each selected tree. The tree chambers were installed in June 2017 in the following order: at the bottom part of the tree stem (~10 cm above the soil) and at 80 and 170 cm above the ground. The rectangular shape stem chambers were made of transparent plastic containers, including removable airtight lids (Lock & Lock Co Ltd, Seoul, Republic of Korea). For the chamber, preparation see Schindler et al.54. Two chambers per profile were set randomly across 180° and interconnected with tubes into one system (total volume of 0.00119 m³) covering 0.0108 m² of stem surface. A pump (model 1410VD, 12 V; Thomas GmbH, Fürstenfeldbruck, Germany) was used to homogenize the gas concentration prior to sampling. Chamber systems remained open between each sampling campaign. During 60 measurement campaigns, four gas samples (each 25 ml) were collected from each chamber system via septum in a 60 min interval: 0/60/120/180 min sequence (sampling time between 12:00 and 16:00) and stored in pre-evacuated (0.3 bar) 12 ml coated gas-tight vials (LabCo International, Ceregidion, UK). The gas samples were analyzed in the laboratory at the University of Tartu within a week using gas chromatography (GC-2014; Shimadzu, Kyoto, Japan) equipped with an electron capture detector for detection of N2O and a flame ionization detector for CH4. The gas samples were injected automatically using Loftfield autosampler (Loftfield Analytics, Göttingen, Germany). For gas-chromatographical settings see Soosaar et al.55.

Soil and stem flux calculation

Fluxes were quantified on a linear approach according to the change of CH4 and N2O concentrations in the chamber headspace over time, using the equation according to Livingston and Hutchison56.

Stem fluxes were quantified on a linear approach according to the change of N2O concentrations in the chamber headspace over time. A data quality control (QC) was applied based on R2 values of linear fit for CO2 measurements. When the R2 value for CO2 efflux was above 0.9, the conditions inside the chamber were applicable, and the calculations for N2O gases were also accepted in spite of their R2 values.

To compare the contribution of soil and stems, the stem fluxes were upscaled to hectares of ground area based on average stem diameter, tree height, stem surface area, tree density, and stand basal area estimated for each period. A cylindric shape of the tree stem was assumed. To estimate average stem emissions per tree, fitted regression curves for different periods were made between the stem emissions and height of the measurements as previously done by Schindler et al.54.

EC instrumentation

EC system was installed on a 21 m height scaffolding tower. Fast 3-D sonic anemometer Gill HS-50 (Gill Instruments Ltd., Lymington, Hampshire, UK) was used to obtain three wind components. CO2 fluxes were measured using the Li-Cor 7200 analyser (Li-Cor Inc., Lincoln, NE, USA). Air was sampled synchronously with the 30 m teflon inlet tube and analyzed by a quantum cascade laser absorption spectrometer (QCLAS) (Aerodyne Research Inc., Billerica, MA, USA) for N2O concentrations. The Aerodyne QCLAS was installed in the heated and ventilated cottage near the tower base. A high-capacity free scroll vacuum pump (Agilent, Santa Clara, CA, USA) guaranteed an airflow rate of 15 L min−1 between the tower and gas analyzer during the measurements. Air was filtered for dust and condense water. All measurements were done at 10 Hz and the gas-analyzer reported concentrations per dry air (dry mixing ratios).

Eddy-covariance flux calculation and data QC

The fluxes of N2O were calculated using the EddyPro software (v.6.0-7.0, Li-Cor) as a covariance of the gas mixing ratio with the vertical wind component over 30-min periods. Despiking of the raw data was performed following Mauder et al.57. Anemometer tilt was corrected with the double-axis rotation. Linear detrending was chosen over block averaging to minimize the influence of possible fluctuations of a gas analyzer. Time lags were detected using covariance maximization in a given time window (5 ± 2 s was chosen based on the tube length and flow rate). While Webb-Pearman-Leuning (WPL) correction58 is typically performed for the closed-path systems, we did not apply it as water correction was already performed by the Aerodyne and the software reported dry mixing ratios. Both low and high-frequency spectral corrections were applied using fully analytic corrections59,60.

Calculated fluxes were filtered out in case they were coming from the half-hour averaging periods with at least one of the following criteria: more than 1000 spikes, half-hourly averaged mixing ratio out of range (300–350 ppb), QC flags higher than 761.

The footprint area was estimated using Kljun et al.62 implemented in TOVI software (Li-Cor Inc.). A footprint allocation tool was implemented to flag the non-forested areas within the 90% cumulative footprint and fluxes appointed to these areas were removed from the further analysis.

Storage fluxes were estimated using concentration measurements from the eddy system (Eq. (1)), assuming the uniform change within the air column under the tower during every 30 min period63,64:

$${\mathrm{S}} = {\Delta}{\mathrm{C}}/{\Delta}{\mathrm{t}} \ast {\mathrm{z}}_{\mathrm{m}},$$
(1)

where S is storage, ΔC is change in the dry mixing ratio of N2O, Δt is time period (30 min), zm is measurement height (21 m).

In the absence of a better estimate or profile measurements, these estimates were used to correct for storage change. Total flux values that were higher than eight times the standard deviation were additionally filtered out (following Wang et al.36). Overall, the QC procedures resulted in 61% data coverage.

While friction velocity (u*) threshold is used to filter eddy fluxes of CO265, visual inspection of the friction velocity influence on N2O fluxes demonstrated no effect. Thus, we decided not to apply it, taking into account that the 1–9 QC flag system already marks the times when the turbulence is not sufficient.

To obtain the continuous time-series and to enable the comparison to chamber estimates over hourly time scales, gap-filling of N2O fluxes was performed using marginal distribution sampling method implemented in ReddyProcWeb online tool (https://www.bgc-jena.mpg.de/bgi/index.php/Services/REddyProcWeb) (described in detail in Wutzler et al.66).

MATLAB (ver. 2018a-b, Mathworks Inc., Natick, MA, USA) was used for all the eddy fluxes data analysis.

Ancillary measurements

Air temperature, relative and absolute humidity were measured within the canopy at 10 m height using the HC2A-S3—Standard Meteo Probe/RS24T (Rotronic AG, Bassersdorf, Switzerland) and Campbell CR100 data logger (Campbell Scientific Inc., Logan, UT, USA). The potential amount of dissolved N2O in the atmospheric water was calculated based on the absolute humidity and the maximum solubility of N2O in water67. DPD was calculated from air temperature and estimated dew point temperature to characterize the chance of fog formation within the canopy. The solar radiation data were obtained from the SMEAR Estonia station located at 2 km from the study site68 using the Delta-T-SPN-1 sunshine pyranometer (Delta-T Devices Ltd., Cambridge, UK). The cloudiness ratio was calculated as the ratio of diffuse solar radiation to total solar radiation.

Near-ground air temperature, soil temperature (Campbell Scientific Inc.), and SWC sensors (ML3 ThetaProbe, Delta-T Devices, Burwell, Cambridge, UK) were installed directly on the ground and 0–10 cm soil depth close to the studied tree spots. During six campaigns from August to November 2017, composite topsoil samples were taken with a soil corer from a depth of 0–10 cm for physical and chemical analysis using standard methods69.

Statistical analysis

R version 4.0.2 (R Development Core Team, 2020) was used to examine, analyze and visualize the data. The significance level (alpha) considered for all the tests was 0.05. The “akima” package version 0.6–2.1 was used to create interpolated contour plots representing a three-dimensional surface70 by plotting soil temperature and SWC against soil N2O emissions as the independent variable. Linear regression models were fitted and Spearman’s rank correlation coefficients were shown for change of SWC and soil N2O flux in period drought onset and air temperature and soil N2O flux in period freeze-thaw. Spearman’s rank correlation coefficients were also shown characterizing the relationship between the monthly average number of days with a high chance of sunshine and fog formation and the difference between the N2O flux from soil and ecosystem. Regarding all measurements of soil temperature, SWC, and soil N2O flux, relationships were better represented by nonlinear than linear models. In addition, the Bragg equation with four parameters71 was used for describing the relationship between SWC and soil N2O flux in the dry period. A workflow for the nonlinear regression analysis was used72 and regression models were fitted in R using functions lm, nls, or loess.