Study area
The study was conducted in post-agricultural forests at different successional stages within and around the Yoko Forest Reserve. The Yoko Forest Reserve (N00°17′; E25°18′; mean elevation 435 m a.s.l.) is located approximately 29 to 39 km southeast of Kisangani in DR Congo (Fig. 1). The vegetation of the region is semi-deciduous rainforest and the climate falls within the AF type in the classification of Köppen. The annual rainfall volume ranges between 1418 mm and 1915 mm and the annual mean temperature is 24.2°C (De Ridder et al. 2014). We collected the samples in a setup of 18, 40 by 40 m plots, which were set up along a secondary succession gradient that consists of six successional stages, i.e., three plots per successional stage, beginning with agricultural fields (AG) and ending with an old growth forest (OG). The AG sites consisted of cassava fields that were prepared (e.g., the forest on the site was cleared and burned) in May 2018 and planted in July 2018. The secondary forests were aged 5, 12, 20 and 60 years-old. Soils in the region are highly weathered Ferrasols, being poor in nutrients, with low pH and dominated by a sandy texture (Table 1).
Table 1
Vegetation structure and soil characteristics of the plots at the study site. Basal area (BA), aboveground biomass (AGB), bulk density (BD) and soil organic carbon (SOC) are, respectively, basal area, aboveground biomass, soil bulk density and soil organic carbon. AG, 5yrs, 12yrs, 20yrs, 60yrs and OG stand for the agricultural site, 5, 12, 20, 60 years-old secondary forest and old growth forest, respectively. Means, standard deviations (sd) and Kruskal-wallis P-values are reported.
|
AG
|
5yrs
|
12yrs
|
20yrs
|
60yrs
|
OG
|
|
|
mean (sd)
|
mean (sd)
|
mean (sd)
|
mean (sd)
|
mean (sd)
|
mean (sd)
|
P-value
|
Stem density (stem ha− 1)
|
-
|
606.3 (39.0)a
|
456.3 (43.8)abc
|
564.6 (140.6)ab
|
410.4 (74.0)bc
|
395.8 (52.0)c
|
0.04
|
Plot BA (m2 ha− 1)
|
-
|
15.9 (6.4)a
|
22.7 (3.4)ab
|
23.7 (4.9)ab
|
25.7 (1.8)ab
|
35.7 (2.9)b
|
0.04
|
Plot AGB (Mg ha− 1)
|
-
|
42.5 (17.4)a
|
79.8 (17)ab
|
96 (28.1)ab
|
182 (10.4)b
|
427.8 (105.8)c
|
0.012
|
Clay (%)
|
6.4 (4.5)
|
8.8 (7.4)
|
9.1 (5.1)
|
8.3 (8.2)
|
12.8 (10.1)
|
14.1 (6.4)
|
0.32
|
Silt (%)
|
6.9 (5.5)
|
7.5 (5.0)
|
8 (5.9)
|
9.3 (9.2)
|
5.9 (6.2)
|
5.3 (3.3)
|
0.66
|
Sand (%)
|
86.7 (7.5)
|
83.7 (10.8)
|
82.9 (6.2)
|
82.4 (11.4)
|
81.3 (10.4)
|
80.5 (6.4)
|
0.49
|
BD (g cm− 3)
|
1.3 (0.1)ab
|
1.1 (0.2)bc
|
1.3 (0.1)a
|
1.2 (0.1)abc
|
1.2 (0.1)abc
|
1.1 (0.1)c
|
0.034
|
SOC (%)
|
0.7 (0.3)
|
1.1 (0.8)
|
0.5 (0.3)
|
0.9 (0.3)
|
0.8 (0.5)
|
3.3 (1.0)
|
0.43
|
pH H20 top soil (0–10 cm)
|
-
|
4.7 (0.1)a
|
4.8 (0.2)a
|
4.7 (0.1)ab
|
4.5 (0.1)c
|
4.5 (0.1)bc
|
0.002
|
Nitrogen deposition and leaching: sampling protocol, sample and data analysis
The bulk precipitation, throughfall and lysimeter samples were collected weekly during one hydrological year (i.e. from July 2018 to August 2019) in the three agricultural plots and the twelve secondary forest plots (i.e. 5yrs, 12yrs, 20yrs and 60yrs). The OG plots were not sampled for N deposition and leaching because of logistical constraints but these plots were already studied earlier (see Bauters et al. 2018, 2019). The bulk precipitation and throughfall samples were collected using polyethylene (PE) funnels, supported by a wooden pole of 1.5 m height that was connected to a 5-L PE container with a PE tube. A nylon mesh was placed in the neck of the funnels to avoid contamination of the collected water by large particles. The containers were buried and covered by litter to avoid the growth of algae and to keep the samples cooler. In total, we used a network of 144 collectors (i.e. 120 collectors for 15 plots of throughfall in the forest plots and 24 collectors of 3 plots of bulk precipitation in the agricultural plots). The setup of bulk and throughfall collectors consisted of a rectangle of 8 by 24 m with 8 m distance between each collector.
Soil solution at 80 cm depth was sampled with 4 lysimeters per plot, installed in the vicinity of the throughfall and bulk precipitation collectors. These ’suction cup’ lysimeters consisted of a PVC tube fitted with a porous ceramic cup (Eijkelkamp Soil and Water, Giesbeek, The Netherlands) and connected to a buried opaque 2-L glass bottle by a PE tube. Using a battery-powered vacuum pump (Prenart Equipment, Copenhagen, Denmark), we applied an underpressure of -500 hPa at each sampling occasion.
For the throughfall, open field rainfall, and lysimeter sampling, the weekly sampling protocol was identical. As such, the water volume of each collector was measured in the field and the materials (i.e. recipients, funnels, mesh) were replaced and subsequently rinsed with distilled water. A plot-level composite sample of each collector category was made, considering the volume of each individual sample. These samples were transported to the lab in Kisangani, filtered with a 0.45 µm nylon membrane filter and put in the freezer for storage within 24 hours. Finally, the samples were shipped in batch to Ghent University in Belgium for chemical analysis. The NH4+ concentration was determined colorimetrically by the salycilate-nitroprusside method (Mulvaney 1996) on an autoanalyzer (AA3; Bran and Luebbe, Norderstedt, Germany) while NO3− was measured by Ion Chromatography (Thermo-Scientific, Pittsburgh, Pennsylvania, USA). Due to some technical problems a bunch of samples was not analyzed for NH4+. Total dissolved nitrogen (TDN) was determined in the water samples with the persulfate oxidation method. For this, an oxidizing solution containing NaOH, H3BO3 and K2S2O8 was added to the sample (Yasui-tamura et al. 2020), which was subsequently placed in an autoclave at 121°C for 1 hour in order to convert NH4+ and dissolved organic N (DON) into NO3−. We calculated the DON concentration as the difference between TDN and the sum of NH4+ and NO3− analyzed in a non-digested sample.
The water flux for bulk precipitation and throughfall was calculated by dividing the average water volume by the surface area of the funnel. We calculated the intercept evaporation at plot scale as the difference between rainfall and throughfall under the assumption that the stem flow was negligible (Zimmermann et al. 2013). Nitrogen deposition fluxes were calculated by multiplying the water volumes with the corresponding concentration of N species. Further, we used the filtering approach of the canopy budget model to estimate dry deposition and canopy exchange from the throughfall and bulk precipitation data (Ulrich 1983; Van Langenhove et al. 2020). Consequently, the total atmospheric deposition was obtained by summing the bulk and the dry deposition, and the canopy exchange by subtracting total atmospheric deposition from throughfall deposition. The N leaching flux at 80 cm depth was calculated by the chloride mass balance (CMB) method (De Schrijver et al. 2004). This method is based on the assumption of mass conservation between the input of atmospheric chloride and the chloride flux in the subsoil (Eriksson and Khunakasem 1969).
Soil N metrics
We sampled topsoil (0–10 cm depth) at five locations and combined them into one composite sample per plot. The samples were oven dried at 60°C for 48h, homogenized and analyzed. The N content and δ15N was determined using an elemental analyzer (automated nitrogen carbon analyzer; ANCA-SL, SerCon, UK), interfaced with an isotope ratios mass spectrometer (IRMS; 20–22, SerCon, UK). The soil texture was determined by the pipet method (Burt, 2004). After removal of the organic matter with H2O2 30%, soil particles were dispersed by Na4P2O7 (aq). The total soil N stock was calculated as the product of soil N content, bulk density, and the soil upper 10 cm, and extrapolated to the plot-scale surface area. We extracted NH4+ and NO3− with 1M KCl. For each sample, 30g soil was extracted with 60 ml 1M KCl, shaken for 1h and filtered (MN615; Macherey-Nagel, Darmstadt, Germany). The NH4+ and NO3− concentrations in the KCl extracts were determined colorimetrically using an auto analyzer (AA3; Bran and Luebbe, Norderstedt, Germany).
Nitrous oxide (N2O) emissions were determined every two weeks during one hydrological year (e.g. from April 2019 to May 2020 because of delayed arrival of static chambers at the study site) using manual static chambers (Hutchinson and Mosier 1981). The chamber bases were installed permanently to avoid soil disturbance. We used chambers made of PVC (diameter = 0.3 m, height = 0.3 m), equipped with an airtight lid, thermocouples, sampling ports, and a vent tube to avoid pressure disturbances. On every sampling occasion, chambers were closed for 1h and the headspace samples taken at 20 minutes interval using a 20 mL syringe and stored in pre-evacuated 12mL containers (Exetainer, Labco, Lampeter, UK). At each sampling time point, temperature inside the chamber was measured using a thermocouple (Type T, Omega Engineering, Stanford, CT, USA). The samples were analyzed at ETH Zurich for N2O using gas chromatography (456-GC; Scion Instruments, Livingston, UK). We calculated the N2O fluxes according to Eq. 1.
F= [(V*P)/(R*S*T)]*ΔC/Δt (Eq. 1)
Where V is the volume of static chamber, P is the pressure, R is the gas constant (0.08206 [L atm K− 1 mol− 1]), S is the area, T is the temperature and ΔC/Δt is the slope of linear regression model representing the rate of concentration change in one hour.
Plant N metrics
A tree species inventory was conducted in the 40 by 40 m plots. The diameter at brest height (DBH) of trees with diameter ≥ 10 cm was measured, marked, and identified at species level. Fully expanded, sun leaves were sampled from the canopy for all species that contributed to 85% of the cumulative plot basal area. The leaves were oven-dried at 70°C for 48h and ground before chemical analyses. Leaf C and N content and leaf δ15N were analyzed using an elemental analyzer (ANCA-SL, SerCon, UK), interfaced with an isotope ratio mass spectrometer (20–22, SerCon, UK). We further calculated the community-weighted mean (CWM) for these measured leaf functional traits, weighed with basal area per species per plot.
Litter traps (dimeter = 0.68 m) were installed in parallel to the throughfall funnels. Consequently, as for throughfall, eight litter traps per plots were set up in a rectangular 8 by 24 m sub-plot and were sampled once a week from July 2018 to August 2019 in twelve plots of secondary forests (5yrs, 12yrs, 20yrs and 60yrs). As we were interested in the fine litterfall, branches of diameter > 2 cm were discarded when present in our traps. Subsequently, litter samples were oven dried at 70°C for 24 hours, weighed and ground before chemical analysis. A subset of litter samples, selected to be representative of seasons, was analyzed for N content using an elemental analyzer (ANCA-SL, SerCon, UK). The litter N flux was calculated by multiplying litter mass by litter N content.
We used ingrowth cores to assess fine roots production from December 2019 to February 2021. The sampling was conducted with 3 month intervals according to an established protocol (Metcalfe et al. 2007; Marthews et al. 2014). We washed the fine roots, oven dried them for 24h at 70°C, and finally weighted and ground them before chemical analysis. The amount of root dry mass collected after three months represents the trimestral root production and subsequently the annual production is the sum of the trimestral production. A subset of root samples, selected to be representative of seasons, was analyzed for N content using an elemental analyzer (ANCA-SL, SerCon, UK). The contribution of fine root production to the plot N was calculated by multiplying the fine root mass production by the corresponding fine root N content.
Statistical analysis
We conducted one-way analysis of variance (ANOVA) to test the variation of the different components of the N cycle and metrics among successional stages. Additionally, we applied the Tuckey HSD post-hoc test. In case the underlying assumptions of ANOVA, mainly homoscedasticity and normality, were not met, we used the Kruskal-Wallis and subsequently its post-hoc Dunn’s test to analyze the data. We performed the statistical analyses using the R programming language 3.6.1 (R Development core 2018).