Long-term nutrient addition increased CH4 emission from a bog through direct and indirect effects

Peatlands are globally significant sources of atmospheric methane (CH4). While several studies have examined the effects of nutrient addition on CH4 dynamics, there are few long-term peatland fertilization experiments, which are needed to understand the aggregated effects of nutrient deposition on ecosystem functioning. We investigated responses of CH4 flux and production to long-term field treatments with three levels of N (1.6–6.4 g m−2 yr−1 as NH4NO3), potassium and phosphorus (PK, 5.0 g P and 6.3 g K m−2 yr−1 as KH2PO4), and NPK in a temperate bog. Methane fluxes were measured in the field from May to August in 2005 and 2015. In 2015 CH4 flux was higher in the NPK treatment with 16 years of 6.4 g N m−2 yr−1 than in the control (50.5 vs. 8.6 mg CH4 m−2 d−1). The increase in CH4 flux was associated with wetter conditions derived from peat subsidence. Incubation of peat samples, with and without short-term PK amendment, showed that potential CH4 production was enhanced in the PK treatments, both from field application and by amending the incubation. We suggest that changes in this bog ecosystem originate from long-term vegetation change, increased decomposition and direct nutrient effects on microbial dynamics.

The emission of CH 4 from peatlands is the balance between production of CH 4 by methanogenic archaea which require anoxic conditions, and oxidation of CH 4 by methanotrophic bacteria mainly in the oxic layer. When nutrient addition leads to increased CH 4 efflux or production, this has often been attributed to the increased abundance and growth of vascular plants, specifically sedges (Cyperaceae), leading to more efficient gas transfer from the rooting zone to the atmosphere and an increased supply of fresh substrates to the methanogenic zone compared to moss-dominated vegetation [21][22][23] . Overall, better nutrient availability alleviates nutrient limitation of microbes and can increase both production and oxidation of CH 4 20-33 . Methane efflux may also increase because of inhibition of CH 4 consuming bacteria by NH 4 . On the other hand, NO 3 can suppress CH 4 production through competition for substrate by denitrification and it can be used as an oxidant increasing CH 4 consumption, thus changing the emission rate 26 . Long-term field experiments may foster our understanding of aggregated effects of nutrient manipulation on ecosystem functioning. There are, however, few long-lasting fertilization experiments in peatlands. In one of these, an oligotrophic fen, increased potential CH 4 production was attributed to vegetation change resulting in increased substrate supply and larger methanogenic population 21 .
We studied CH 4 dynamics at long-term fertilization experiments at Mer Bleue bog, Ontario, Canada, where treatments with control and three levels of N as NH 4 NO 3 with or without constant PK (KH 2 PO 4 ) addition started in years 2000 (experiment I with N, PK, NPK) and 2005 (experiment II with N only) (Fig. 1). We measured CH 4 fluxes between the ecosystem and the atmosphere weekly from May to August in 2005 (treatment years 1 and 6, for experiments I and II, respectively) and 2015 (treatment years 11 and 16). We incubated peat samples collected from the fertilization plots in 2016 to examine variations in the rates of potential CH 4 production and consumption and supplemented this by examining the effect of laboratory PK addition on potential CH 4 production. To quantify the peat subsidence visible through relative water table, we use spatially representative GPS data on peat surface elevation collected in years 2011 and 2013.
We hypothesized that 1) treatments with the highest long-term NPK loading would show the largest CH 4 emission driven by a change from Sphagnum moss to dwarf shrub dominance associated with a faster litter decomposition rate and a rise in the WT; 2) potential CH 4 production rates in peat samples from the field treatments would be higher in nutrient treatments with largest vegetation changes owing to reduced microbial nutrient limitation; and 3) potential CH 4 production rates would increase with a laboratory PK amendment, with a larger response in the N-only field treatments than in NPK treatments, because of their nutrient limitation relative to the increased N supply.

Study site and experimental set-up
The study was conducted in a long-term fertilization experiment at Mer Bleue bog, an ombrotrophic peatland complex covering 28 km 2 in eastern Ontario, Canada (45.410017 N, 75.518348 W) 34 . The mean annual temperature (1971-2000) was 6.0 °C and the mean annual precipitation was 943 mm, of which 351 mm fell between May and August 35 . Estimated wet atmospheric deposition of N ranges from 0.5 to 0.8 g N m −2 yr −1 based on Canadian Time frame and major responses in vegetation, ecosystem CO 2 exchange (net exchange, photosynthesis, and respiration), and water table relative to control in the fertilization experiments I and II at Mer Bleue bog. Vertical arrows indicate the years when a response was documented. *Treatments 3.2NPK and 6.4NPK were started in 2001, 1 Bubier et al. 41 , 2 Juutinen et al. 5 , 3 Larmola et al. 7 . and US deposition observation networks (see also Turunen et al. 36 ) and is small compared to the global N deposition gradient 8 .
The fertilization experiment is located in the northwestern part of the bog, where the surface has only a slight (20 to 30 cm) hummock and lawn topographic variation with hummocks comprising about 70% of the area 37 38,39 . New species that have appeared in the long-term fertilized plots include the fern Thelypteris palustris (L.) A.Gray, and numerous moss species, for example Pleurozium schreberi (Brid.) Mitt. and Aulacomnium palustre (Hedw.) Schwägr 7 .
The fertilization experiment includes two sets of treatments 5,7,[40][41][42] . The first experiment (I) was established in 2000-2001, which includes a control (C1), low level N-only, PK, and NPK additions at three N levels (1.6, 3.2 and 6.4 g N m −2 yr −1 as NH 4 NO 3 ) (Fig. 1, Table 1). These levels range from about double the ambient N deposition at Mer Bleue to the high deposition rates occurring, for example, in Europe 8 . The second experiment (II) was started in 2005 to include a control (C2) and the 3.2 and 6.4 g N m −2 yr −1 levels without PK. PK treatment was applied at levels 0 and 5.0 g P and 6.3 g K m −2 yr −1 as KH 2 PO 4 ( Table 1). Each treatment has three plots, which are 3 m × 3 m and have at least 1 m buffer zone between them. The nutrients were added dissolved in distilled water corresponding to precipitation of 2 mm as seven applications every three weeks between May 1 st and August 31 st . The controls received the same amount of distilled water.  Fig. 2b). Following the typical seasonal pattern, the CH 4 fluxes at the control plots (C1) rose from May onwards as the peat warmed and peaked in late July in both years, with a slightly larger seasonal peak in 2005 than in 2015 (Fig. 2c).

Results
Detrended correspondence analysis (DCA) revealed changes in plant community composition between treatments and between years 2005 and 2015 (Fig. 3). Plots in treatments 3.2NPK and 6.4 NPK diverged from the other plots in 2005 along the first axis, and this variation in plant community composition became larger by 2015. In addition, plots in treatments 1.6NPK and 6.4N moved along the first and second axes over the 10 year period. The first DCA axis represented gradients in cumulative N loading, PK addition and CH 4 flux based on a post-hoc fit of environmental variables in species space. The second axis represented gradient in WTD and year. Methane flux correlated with cumulative N (r = 0.531), PK (r = 0.426) and WTD (r = −0.355), cumulative N loading correlated also with WTD (r = −0.557).
In 2005, the experiments had been exposed to treatments for 1 and 6 yr and Sphagnum moss cover decreased in all six year old nutrient addition treatments. In 2015, Sphagnum also decreased in the nutrient treatments started in 2005 (Fig. 3b,c), and abundances of both Sphagnum and Polytrichum strictum mosses correlated negatively with cumulative N load (Fig. 3b). An increase in vascular plants, mainly dwarf shrub abundance, was synchronous with a decrease in Sphagnum coverage, but with a smaller magnitude (Fig. 3b,c). Differences in WTD between the treatments were apparent only in 2015 when all NPK treatments (16 yr) and 6.4 N treatment (11 yr) had a WTD closer to the peat surface than the rest of the treatments (Fig. 3d,e).
The surface elevation showed effectively no change between 2011 and 2013 for the control plots where the average (±1 SD) difference was 1 ± 3 cm for C1 and 0 ± 2 cm for C2 (Supplementary information Table S1). The peat surface decreased, on average, by 5 ± 3 cm in the 6.4NPK plots and 6 ± 4 cm in the 6.4 N plots. In 2013, the  average surface of the 6.4NPK plots was 13 cm lower than in C1 plots, and that of 6.4 N plots 5 cm below the surface of C2 plots.
Treatment effects on CH 4 flux. Individual CH 4 flux values were averaged for each plot over the study period (May-August) to compare treatments and control. In 2005, the mean CH 4 flux was largest in the PK treatment (22.5 mg m −2 d −1 ), but with no significant treatment effects. In 2015, CH 4 flux in the 6.4NPK treatment was large and differed significantly from the C1 (50.5 vs. 8.6 mg m −2 d −1 , Tables 2 and 3, Fig. 4). In 2005, the relationship between CH 4 flux and WTD was not significant, in part because of the small range in WTD: 22 cm among the 24 collars and 9 cm among the 8 treatments means (Fig. 4). In 2015, the relationship was significant with a wider range in WTD (42 cm among the 27 collars and 24 cm among the 9 treatments means) and the relationship was also significant when data from 2005 and 2015 were combined.
The effect of long-term field treatments and short-term PK amendment on CH 4 production potential. Sampling of peat relative to WTD in 2015 revealed that anaerobic CH 4 production potentials in    Fig. 1). In 2016, peat samples collected from the C1, PK, 6.4 N and 6.4NPK field treatments at a depth of average WT position, had anaerobic CH 4 production potentials ranging from 2 to 22 µg g −1 d −1 (Fig. 5a), with significantly larger rates in the long-term PK field treatments than in C1 and 6.4 N field treatments ( Table 4, Fig. 5b). The laboratory PK amendment showed an increase in anaerobic CH 4 production potentials across the field-treatments (p = 0.086, Table 4, Fig. 5c), with increases, relative to unamended samples, of 16%, 367%, 186% and 256% in the field treatments 6.4 N, C1, 6.4 NPK, and PK, respectively. The degree of decomposition of peat samples was higher by one unit in the Von Post scale in the peat cores of 3.2NPK, 6.4 NPK and 6.4N treatments relative to the controls.

Discussion
In agreement with our hypothesis, CH 4 flux was nearly 5 times greater in the 6.4NPK treatment than in the control in the 16 th treatment year (2015), whereas differences were not significant in 2005. Conditions resulting in increased CH 4 fluxes developed over the long-term and were detected only in the treatment with the longest and largest N loading with PK. The trends in vegetation, WTD, and peat subsidence in other NPK and 6.4N treatments, however, may indicate a similar response in the future. We interpret the response of CH 4 flux to be the result of multiple environmental changes, initiated by a drastic decrease in Sphagnum moss abundance and concurrent increase in dwarf shrub biomass altering the quantity and quality of labile carbon input (Figs 1 and 3). These changes have led to increased ecosystem and heterotrophic respiration indicating stimulated decomposition 7,43 , and to peat subsidence and rise in the WT relative to the peat surface in the NPK plots (Fig. 4). The largest CH 4 fluxes were associated with the highest WT positions, following the general pattern observed among and within peatlands 18 , but both the largest CH 4 fluxes and highest WT position were associated with the largest additions of NPK. Guo 44 found a decrease in methanogen population size and no change in methanotroph populations in response to nutrient treatments at the depth of 30-40 cm in peat samples from the fertilization experiment collected in 2013. There were, however, broad changes in other microbial groups, for example a decrease in the fungi to bacteria ratio, and increased bacterial numbers and diversity in response to nutrient addition. The microbial results were similar to vegetation responses, in that NPK treatments induced stronger responses than N alone due to NP co-limitation at the site 45 Table 3. Analysis of variance on the effects of treatments on log 10   production potential as shown by Robroek et al. 46 . The results of their vegetation manipulation study, contrary to ours, showed that an intact Sphagnum layer with graminoid and ericoid removal resulted in a decrease in CH 4 production potential by reducing root derived substrates. In our study, both the long-term field and short-term lab addition of PK were associated with larger CH 4 production potential than control and N-only treatments (Table 4, Fig. 5). This may be due to better availability of suitable substrates driven by increased plant productivity, changes in litter quality, and stimulation of decomposition 7,43,46 . The immediate increase of CH 4 production potential in response to PK amendments suggests a direct positive nutrient impact on the microbial community, relieving methanogenic and syntrophic microbes from nutrient deficiency. This result agrees with Pinsonneault et al. 43 who found increased activity of β-D-glucosidase enzyme with NPK fertilization at Mer Bleue bog. We cannot distinguish whether methanotrophs and CH 4 consumption were inhibited by higher concentrations of NH 4 owing to its similarity with CH 4 and inhibition of the enzyme responsible for CH 4 oxidation 20,27,31 .
Previous data on CH 4 production from fertilization experiments varying in site characteristics and experimental length show mixed results, illustrating that methanogenic responses to nutrients are complex. Keller et al. 27 found in short-term and 6 yr-long treatments that P and NP (NH 4 , up to 10 g N m −2 yr −1 and 2 g P m −2 yr −1 ) stimulated CH 4 production in a rich fen, but inhibited it in a nutrient poor bog, suggesting that initial conditions, including nutrient availability, may determine the response. Eriksson et al. 21 inferred from long-term data that the CH 4 responses to nutrients are dependent on the WTD and vegetation, and are related to the amount, development and nature of microbial populations. Site properties, thus, may determine whether nitrate and other N oxides decrease CH 4 production in N-only treatments if methanogenesis becomes outcompeted by energetically more efficient terminal electron acceptors and competition for substrate by denitrification 26,47,48 .
Mer Bleue, a bog dominated by dwarf shrubs and mosses with a deep WT, has not exhibited a shift to graminoid vegetation associated with an increase in CH 4 fluxes as in some other N addition experiments. In a boreal ombrotrophic bog, 6 years of NH 4 NO 3 addition (10 g N m −2 yr −1 ) resulted in a doubling of CH 4 emission,  Table 4), open symbol indicate treatment with distilled water (DW) and closed with PK solution (PK). Significant differences between treatments (Bonferroni multiple comparisons) are indicated in the b panel.  Table 4. Three-way nested analysis of variance on the effects of treatment, plot nested within the treatment, PK amendment and the interaction treatment × PK amendment on potential CH 4 production rates. attributed to the increase in sedge Eriophorum vaginatum abundance 23 . In Degerö Stormyr, an oligotrophic fen in Sweden, N addition decreased CH 4 fluxes in plots with a high sedge cover, but increased fluxes slightly in plots with a small sedge cover; this was interpreted as being caused by changes in the distribution of roots and consequent effects on substrate supply 49 . Later in the same experiment, 12 years of fertilization with 3 g N m −2 yr −1 increased potential CH 4 production significantly in the layer of maximal CH 4 production near the water table. This was explained by a marked increase in the abundance of sedge vegetation and increases in the amount of labile substrate from root exudation 21 . Similar to our treatments with N-only, no or very small effects have been found in short-to mid-term NH 4 or NH 4 NO 3 addition experiments in oligotrophic fens 22,24,50 .
In conclusion, elevated CH 4 fluxes at Mer Bleue are the result of several ecosystem processes altered by 16 yr of fertilization at a rate of 6.4 g N m −2 yr −1 , plus P and K. The long-term addition of NPK has affected microbial community structure 44 and enzyme activities 43 indirectly through changes in vegetation and peat properties. Wetter conditions resulting from loss of moss, increased decomposition and peat subsidence in the highest NPK treatment have moved the zone of CH 4 production closer to the peat surface, thus increasing the availability of added nutrients and substrates to methanogens, and likely reducing the proportion of CH 4 consumption. In addition, PK addition seems to stimulate CH 4 production directly, while the addition of NH 4 NO 3 alone has had little direct effect on either CH 4 production or consumption. However, comparison of this study with the literature shows that responses of CH 4 dynamics are not universal across peatland ecosystems and thus merits further long-term investigation of multiple processes in the ecosystem.

Methods
Methane flux measurements. Methane fluxes were measured approximately weekly from the end of May to the end of August in 2005 and 2015 using a non-steady state, static, closed chamber method employing permanent collars with an area of 0.045 m 2 inserted into the peat to a depth of 15-20 cm in each of the treatment plots. We used opaque chambers made of polycarbonate bottles (volume 18 L), which were covered with aluminum foil to reduce heating. We randomized the order in which the plots were measured between 8 am and 4 pm. In total 15 and 24 measurements per collar were made in 2005 and 2015, respectively; regardless of these differences, they were equally distributed over the May through August measurement period and the means should not be biased.
The top of the collars were grooved and were carefully cleaned of litter and filled with water to ensure airtight conditions before placing the chamber. Chamber headspace air was sampled 6 times at 0, 5, 10, 15, 20, and 30 min. Air samples (30 mL) were withdrawn through tubing and a three-way stopcock after mixing the air in the chamber and tubes by pumping with a 60 mL plastic syringe. The CH 4 concentration was determined within a day after collection in 2005. In 2015, the air sample was injected into a 12 mL pre-evacuated glass vial and analyzed within 1-2 weeks. Methane concentrations were analyzed using a Shimadzu Mini II gas chromatograph (Shimadzu, Kyoto, Japan) and Varian CP-3800 gas chromatograph (Varian, Inc., Walnut Creek, CA, USA) in 2005 and 2015, respectively, both equipped with a packed column (Poropak Q column with 80/100 mesh, Alltech, 8 Deerfield IL, USA), flame ionization detector (FID) SRI-310C, with column and detector temperature of 50 °C and 110 °C, respectively. The gas chromatographs were calibrated for CH 4 with multiple standards (0-22 ppm).
Fluxes of CH 4 were calculated from the slope of the linear regression of gas concentration in the chamber headspace against time, chamber volume corrected for collar volume, surface area, and chamber temperature and ambient pressure. Each measurement series was checked visually by plotting CH 4 concentration against time to check for leakage, ebullition or saturating CH 4 concentration and unrealistically large initial CH 4 concentration. Less than 5% of measurements were discarded and the measurements retained for analysis were associated with time series regression coefficients of determination of 0.90 or greater if there was an increase or decrease in the concentration. This rule was not followed in the case of near zero fluxes, i.e. when there was no change in the concentration over the time.
Characterization of surface elevation, water table position, temperature and vegetation. We had observed the 6.4 NPK plots to be wetter than surrounding plots and to have subsiding peat surfaces 5,7 . To quantify these changes, a systematic survey of surface elevation was conducted in a subset of treatments (C1, C2, 6.4NPK and 6.4N) in 2011 and 2013. Calibration points as well as precise location points and their respective elevation above sea level were measured in treatment plots using a differential global positioning systems (DGPS) GRX-1 (Sokkia Corp, Mississauga, Ontario) and Trimble 5800 (Trimble Navigation Ltd, Dayton, Ohio). The plots were surveyed from a movable ladder set between the boardwalk and a fence with no disturbance to the peat and vegetation. The DGPS was attached to a levelled rod with a piece of plywood to gently keep it at the moss surface for measuring the elevation of the peat or Sphagnum moss surface. The systems were referenced to a base station that was positioned within 500 m from the surveying area on stable ground or over a Natural Resources Canada monument and permanent benchmark points near the plots were used to calibrate between the two years. Location and elevation above the sea level were recorded on approximately 220 grid points per plot in 2011. In 2013, approximately 104 grid points were measured in C1 and 6.4 NPK plots and 33 grid points were measured in the less spatially variable C2 and 6.4N plots. The horizontal and vertical instrumental precision were ±0.01 m and ±0.015 m (GRX-1, Sokkia) and ±0.01 m and ±0.02 m (Trimble 5800), respectively. ArcMap v10.1 (ESRI, 2013) was used to create surface rasters for the 2011 and 2013 years. The difference between raster surfaces for these two years was used to identify the change in elevation within the plots.
Water table depth (WTD) and chamber temperature were measured concurrent with CH 4 flux measurements. Water table position was measured by determining its distance from the peat surface in perforated 2 cm diameter tubes that were inserted next to each collar. To characterize the year-to-year differences in conditions affecting CH 4 fluxes, we also used continuous air and peat temperature and WTD measurements, averaged to daily values, made at the eddy covariance flux tower 100 m south of the fertilization plots. The abundance of vegetation and species composition was characterized in 60 × 60 cm quadrats inside the treatment plots by visually estimating the SCiEntiFiC REPORTS | (2018) 8:3838 | DOI:10.1038/s41598-018-22210-2 % cover of species (2005) and using the point intercept method (2015) 51 . The point intercept values were scaled to % values. We present here data measured in July, representing approximately peak above-ground biomass.
Potential CH 4 production rates. We conducted two laboratory incubations of peat collected from the fertilization plots to identify 1) the effect of sample location relative to WT on potential CH 4 anaerobic production and aerobic consumption rates across all field treatments, and 2) the effect of selected field treatments (C1, C2, PK, 6.4 N, and 6.4NPK) and laboratory amendment with PK on potential anaerobic CH 4 production. The first incubation experiment was a pre-study to help focus sampling for the second incubation experiment and, as it showed a potential PK effect, we conceived the second incubation experiment. For the first incubations one 40 cm long 10 cm × 10 cm square peat core was collected from each plot in May 2015 and split into 10 cm sections. The depth from the surface and the distance of each section to the WT at the time of sampling was recorded (at the main flux tower weather station, water table was about 33 cm below the hummock surface). The samples were stored at 5 °C prior to incubation that took place in June 2015.
Approximately 4 to 5 samples per treatment, totaling 38 and 48 samples, were selected for anaerobic and aerobic incubations, respectively. Three subsamples of ~12.5 g of wet peat from the each 10 cm sections were placed in 125 mL Erlenmeyer flasks for both the anaerobic and aerobic incubations and stored in the dark at 20 °C 52 . The anaerobic set was purged with Ultra High Purity nitrogen (N 2 ) for 1 h to ensure anoxic conditions, capped with a rubber SubaSeal stopper and incubated for 21 d. Flask headspace was sampled (3 mL) immediately after purging and every third day, the removed gas being replaced with N 2 . For aerobic CH 4 consumption potential, 125 μL of pure CH 4 was added to ambient air in the flask. Flask headspace was sampled with a 3 mL syringe initially and after 1, 3, and 5 d, being replaced by ambient air. Methane concentrations were determined within a few hours of sampling on a Shimadzu Mini II gas chromatograph as described above. Rates of potential CH 4 production and consumption were calculated from changes in headspace CH 4 concentration, accounting for changes associated with dilution over time, and expressed per mass of dry peat, determined by oven drying peat samples at 50 °C at the end of the incubation. The first incubation study showed that greatest CH 4 production and consumption potential rates were from peat sections just below and at the WT, respectively (Supplementary information Figure S1).
For the second incubation study we sampled the C1, PK, 6.4 N, and 6.4NPK treatments in July 2016 when the WT was 20-60 cm below the surface of the fertilization plots. One sample (10 × 10 × 10 cm) from each treatment plot was collected at the WT, the zone of maximum CH 4 production according to pre-study, and divided into 10 subsamples of 12.5 g for anaerobic incubation (prepared as described above). One subsample from each plot was oven-dried to estimate the moisture content. Five subsamples were amended with PK, by adding 100 μg of P per gram of dry soil, dissolved in water as 0.439 g L −1 KH 2 PO 4 and another five subsamples were treated with deionized water as a control. The flasks were sampled for headspace CH 4 concentration over a 12 d incubation (as the results from the pre-study indicated no change in trend for CH 4 production after 10-12 d), oven dried afterwards and production rate calculated as described above.

Data analyses.
To characterize the relationships between mean CH 4 flux, plant community composition, cumulative N loading and WTD among the fertilization plots, we performed detrended correspondence analysis (DCA) with post-hoc fit of year and environmental variables as supplementary variables. A common analysis for both years was performed by transforming the 2015 vegetation data into percentage values. The DCA was performed on log 10 -transformed centered plant species data using Canoco 5 53 .
Treatment effects on surface topography were assessed by comparing the surface rasters of elevation between years 2011 and 2013 in treatments 6.4 NPK, 6.4 N, C1 and C2 (n = 3 per treatment). Seasonal means of daily CH 4 fluxes and WT were calculated for each plot (n = 3 per treatment, two collars in one of the 6.4NPK plots were first averaged for each day), and used in regression analysis to examine the relationship between WT and CH 4 flux. Log 10 -transformed data were used to fit linear regression allowing a comparison with that for overall Mer Bleue data from Moore et al. 54 . Treatment effects on CH 4 fluxes were analyzed using analysis of variance, handling the two experiments separately. First, treatment and year were used as fixed factors and WT as covariate, and second, when there was significant year × treatment interaction the data were analyzed by years. Differences between the control and treatments were assessed by Dunnett's t post hoc test. Data were log 10 transformed to meet the requirement of equal variances and normal distribution.
Effects of long-term fertilization treatments (control, 6.4N, PK, and 6.4NPK) and laboratory amendment of PK on potential production of CH 4 were analyzed using three-way nested ANOVA with Treatment, Plot nested within the Treatment, and PK amendment as the three factors and the Treatment × PK amendment interaction was also included. Treatment differences were assessed using pairwise Bonferroni comparisons. Five laboratory replicates for each plot and treatment were averaged prior to the analyses (treatment n = 3). Data Availability. The data set generated during the current study will be available through Mount Holyoke College Institutional Archive.