Shifts in soil phosphorus fractions during seasonal transitions in a riparian floodplain wetland

Coppolino, Jacob; Munford, Kimber E.; Macrae, Merrin; Glasauer, Susan

doi:10.3389/fenvs.2022.983129

ORIGINAL RESEARCH article

Front. Environ. Sci., 08 December 2022

Sec. Biogeochemical Dynamics

Volume 10 - 2022 | https://doi.org/10.3389/fenvs.2022.983129

This article is part of the Research Topic The Cold Regions in Transition: Impacts on Soil and Groundwater Biogeochemistry View all 8 articles

Shifts in soil phosphorus fractions during seasonal transitions in a riparian floodplain wetland

Jacob Coppolino¹

Kimber E. Munford¹

Merrin Macrae²

Susan Glasauer¹*

¹School of Environmental Sciences, University of Guelph, Guelph, ON, Canada
²Department of Geography and Environmental Management, University of Waterloo, Waterloo, ON, Canada

Losses of phosphorus from soil to surface waters in agricultural areas have been linked to substantial declines in water quality. Riparian wetlands can potentially intercept phosphorus mobilized from upland soils before it reaches connecting waterways, but the capacity of wetlands to buffer against downstream losses of P is poorly understood, especially in northern temperate zones. In these regions, the spring freshet releases large volumes of water from snowmelt and soil pore water during the time when microbial productivity, which transfers available P into biomass, is low. In addition, losses of P in runoff may be exacerbated by freeze-thaw cycling (FTC) in soil during late winter and early spring through the physical degradation of organic matter. We investigated P dynamics from late fall through spring thaw and into summer to assess P transfers between inorganic, organic and microbial biomass pools, as functions of season and distance from a river. The site is located on the Grand River in southern Ontario, which discharges to Lake Erie, and consists of riparian wetland and wooded areas. Reactive P (Olsen P) and microbial biomass P (P_MBIO) increased with distance from the river and varied more over time in the wetland soil compared to the adjacent wooded area, reflecting higher variability in vegetation, topography and hydrology. The positive correlation between microbial biomass P and microbes linked to ammonification supports the release of N and P through mineralization pathways as spring progresses, with microbial biomass decreasing in June as plant growth increases. There was evidence for leaching of Fe and Al, and lower concentrations of total P, in the transect proximate to the river. Seasonal flooding during spring thaw contributed to a pulse of dissolved reactive P, but temperature monitoring showed that the wetland soil did not experience freeze-thaw cycling. Investigation of FTC using wetland soil in mesocosms indicated that multiple FTC (>3) were necessary to increase the pool of reactive soil P, with the highest amount of soil reactive P observed after six FTC, when dissolved reactive P also tended to increase.

Introduction

Changes in the magnitude and timing of seasonal climatic events are shifting the export of nutrients from farmland into aquatic environments, such as that observed in the Lake Erie watershed (Sims and Sharpley, 2005; Blackwell et al., 2010; Wolf et al., 2017). Excess phosphorus (P) is of particular concern in Lake Erie, the shallowest of the Great Lakes, as it is surrounded by extensive agricultural areas (Wilson et al., 2019). Despite the adoption of best management practices throughout the watershed, harmful and nuisance algal blooms continue (Smith et al., 2019). There is uncertainty on how the changing climate anticipated for the region will impact the mobilization of P from agricultural regions, particularly during winter (Kieta et al., 2018). Thus, an improved understanding of P dynamics is needed in cold agricultural regions. Riparian wetlands are a common Best Management Practice (BMP) used to improve water quality in agricultural areas by trapping nutrients through physical and chemical processes before nutrients reach surface waters (Lowrance et al., 2000; Hickey and Doran, 2004; Rao et al., 2009; Vidon et al., 2019; Young et al., 2021). Among the different types of riparian buffer zones, wetlands can be particularly effective at sequestering P because of their diverse chemical and microbial processes (reviewed in Hoffman et al., 2009). Widespread adoption of any type of BMP by farmers to manage P requires, at a minimum, the proof of efficacy (Wilson et al., 2019), which is currently lacking for riparian wetlands that receive runoff water. One potential concern is that riparian wetland areas are prone to flooding. Although floodplain soils have been shown to retain agricultural P during flooding in warm temperate climates (Gillespie et al., 2018), this has not been confirmed for northern climates that experience substantial spring thaws and consequently have a high potential for the release of dissolved phosphorus from upland soils to surface waters over short periods of time (Liu et al., 2014).

Increased seasonal transfers of P to surface waters during late winter and early spring have been predicted for cold temperate climates due to warmer temperatures and reduced snowpack, which affect freeze-thaw cycling (FTC) (Henry, 2008). FTC in cold temperate soil typically occurs in early spring; it reduces aggregate stability, fosters cell lysis and helps degrade organic matter (Blackwell et al., 2010). The early spring is an especially sensitive time of year for loss of P trapped in pore water because of high runoff, high water levels in rivers, and low biological activity that can trap dissolved P (Blackwell et al., 2010; Kieta et al., 2018). Freeze-thaw cycling disrupts aggregate stability (Wang et al., 2007), which has the potential to affect P mobilization in wetland soil. Such disruption could reduce P mobilization by fostering the sorption of phosphate to organic matter and clay minerals through the exposure of additional sorption sites in response to the mechanical disintegration of soil; alternatively, disaggregation could expose previously occluded P and create new flow paths. As the soil warms and biological activity increases, the bioavailability of phosphorus decreases as P is sequestered by microbes and plants in biomass (Edwards et al., 2006; Guicharnaud et al., 2010; Buckeridge et al., 2016; Xu et al., 2018; Miao et al., 2020). An improved understanding of the impact of FTC on P dynamics is needed to understand the potential for the mobilization of P between riparian wetlands and adjacent surface water bodies during the critical spring-summer period, when a significant proportion of annual P loss occurs (Macrae et al., 2007).

The most reactive or available fractions of soil phosphorus include P associated with orthophosphate, simple polyphosphates, and labile forms of organic phosphate, including microbial biomass P (Dunne and Reddy, 2005). Microbial biomass P is a distinct fraction that can comprise up to 76% of bicarbonate extractable P (Blackwell et al., 2010). Typically, only a small portion of total phosphorus in wetland soil is biologically available (Dunne and Reddy, 2005; Cui et al., 2019). There have been relatively few investigations of microbial P cycling through microbial biomass P in natural wetland soil (Blackwell et al., 2010; Wang et al., 2022). This fraction is of particular interest because of its availability for mineralization (Blackwell et al., 2010). There are even fewer studies that relate microbial community composition to P chemistry in wetland soil, with one study finding no unequivocal correlation between genes involved in P cycling, and P and N fractions (Li et al., 2022). There have been no studies of the spatial distribution of reactive P pools and their response to seasonal spring conditions for natural riparian wetlands in cold temperate agricultural areas. Given the large seasonal fluctuations in water levels for rivers in northern zones, it is important to assess whether riparian wetlands serve as sources or sinks of P, as investigated for warmer climate zones (Noe et al., 2013; Gillespie et al., 2018).

We investigated changes in the fractions of P and in microbial communities across the surface (upper 10 cm) of riparian wetland soil. The wetland that we studied is connected to the Grand River through runoff from surrounding farmland and through spring flooding; this river travels south through southern Ontario to Lake Erie. The objectives were to:

1. Quantify intra-annual changes in reactive P, soil microbial biomass P and total soil phosphorus fractions. We focused on the period from late winter to summer to identify changes in the pools of P during the time of highest biological activity;

2. Relate soil P fractions to spatial location, soil physicochemical characteristics and soil microbial communities;

3. Explore the impacts of freeze-thaw cycling on reactive soil P fractions and microbial activity in controlled experiments to identify how repetitive FTC may shift pools of soil P and impact microbial responses to P.

In the field investigations, we expected to observe shifts in the reactive pools of P especially from early spring to summer that correlate with increased microbial and plant activity. The concentrations of total soil P were anticipated to reflect the diverse moisture and chemical conditions across the site, with drier sites retaining more P. Because we did not observe soil freezing in the field, despite sustained subzero air temperatures, we investigated the impacts of freeze-thaw cycling on soil P using ex situ soil mesocosms.

Materials and methods

Site description

The study area is located along the Grand River in Grand Valley, Ontario [43°58′33.71″N; 80°21′40.48″W] (Figure 1). The wetland is situated in a depression at the base of two gently sloping fields, which are cultivated for hay and fertilized every 5 years with manure. The plant species include characteristic plants for wet soil: red osier dogwood, ostrich fern, spotted joe-pye weed, fragrant bedstraw, swamp milkweed, marsh marigold, Canada anemone, and Canada nettle (Newmaster et al., 1997). Iron mottling was visible in the soil during the spring transition, with nodules around 0.5 cm or smaller. During the summer, the mottles were no longer observed and the soil appeared gleyed below 5 cm. It is a riparian floodplain marsh according to the Canadian Wetland Classification System (National Wetlands Working Group, 2018).

FIGURE 1

FIGURE 1. (A) Location of the site within Ontario (red), (B) Detailed view of site; arrows indicate the flow of water in the Grand River and (C) The sampling grid for wetland and treed zones of the site. Map generated using QGIS. Basemap source: Google Earth.

The wetland is bordered at the northern and southern margins by areas planted with deciduous trees. The area vegetated with trees at the southeast margin of the wetland was included in our study to compare P cycling between the wetland and a more elevated site with distinct vegetation.

Surface runoff at the site travels via drainage ditches between fields, and between fields and the riparian zone. When flows are high, runoff moves through the wetland and treed areas, as well as via flow in the creek between the two zones (Figure 1C). The creek receives additional drainage from fields situated to the northeast (Figure 1). The wetland is covered with snow and ice in the winter and flooded each spring by the Grand River during freshet (Supplementary Figure S1); the treed area is not flooded.

Field sampling of soil and water samples

Soil and water samples were collected on six dates (2020–12–10; 2021–03–12; 2021–03–26; 2021–04–15; 2021–05–05; 2021–06–10). The soil was sampled in three transects (Figure 1) parallel to the Grand River. Ten scoops of soil (around 250 ml each) were collected within an area of 0.5 m diameter at each site and composited. The soil in the treed zone was sampled in a transect at three sites perpendicular to the river and parallel to the wetland sites. Samples were stored at 4°C or air dried immediately. Samples for P_AD and [Metal]_AD (acid digestion) were dried at 105°C. Dried soil was lightly crushed and sieved, depending on the analysis.

Water samples to assess dissolved reactive phosphorus (DRP) were collected at the bank of the Grand River near soil sampling site 2, at the creek outlet, and 50 m further upstream in the creek (Figure 1). An acid washed 1 L Nalgene bottle was filled with water from each of the three sites. Samples were kept at 4°C in the dark and filtered (cellulose acetate, 0.45 µm) within 6 h of collection. Four 10 ml subsamples of the filtrate were analyzed using a Hach portable colorimeter with PhosVer reagent packets (Hach Company). PhosVer reagent packets have a detection range of 0.02–2.50 mg/L DRP. Temperature, moisture, and electrical conductivity were monitored in situ using TEROS 11/12 solar-powered soil probes (METER Group, Oakville, Ontario). Two stations were set up in the wetland for continuous monitoring using a ZL6 data logger with Zentra cloud capabilities, with data recorded every 15 min. The temperature, moisture, and electrical conductivity each had an accuracy rating of; ± 1 from −40°C to 0°C and ±0.5 from 0°C to +60°C, ±0.03 m³/m³, and ±(5% + 0.01 dS/m) from 0 to 10 dS/m ±8% and from 10 to 20 dS/m respectively. Each metric had detection ranges of −40°C to +60°C, 0.00–0.70 m³/m³, and 0–20 dS/m respectively. Each data logger was connected to an ATMOS-14 atmospheric sensor and three TEROS 11/12 soil probes each. The ATMOS-14 sensor monitored air temperature directly above the soil. The ATMOS-14 sensor was rated for accuracy of ±0.2°C and a detection range of –40.0°C to +80.0°C. The TEROS 11/12 soil probes were positioned ∼120° apart at a depth of 10 cm and a distance of ∼5 m from the base of each weather station mast. No soil probes or weather sensors were installed in the treed area because there was insufficient sunlight to power them.

Soil analysis

Air-dried soil was analyzed in triplicate for organic C, inorganic C, and total N by combustion (LECO CN828) at the Food and Agricultural Laboratory, University of Guelph. Soil organic carbon was ∼ 5% in the wetland, with a slightly higher value in early spring and a lower value in late fall. The C:N ratio was ∼11, consistent with the designation of hydric mineral soil (Graf-Rosenfellner et al., 2016; National Wetlands Working Group, 2018).

Phosphorus and other elements in soil were assessed to identify distinct chemical fractions. The Olsen P method uses bicarbonate extraction of soil to determine molybdate-reactive P (P_OI) and “total” inorganic and organic reactive P (P_OT) (Table 1), which includes organic reactive phosphorus and polyphosphates (Olsen et al., 1954; Richardson and Reddy, 2013). To determine Olsen P, 5 g of air-dried and sieved soil (<2 mm) was extracted with 0.5 M NaHCO₃, 0.4 ml of activated carbon (300 g of phosphate-free charcoal with 900 ml of deionized water), and shaken for 16 h (Zhang, and Kovar, 2000). After extraction, the solution was separated by centrifugation and filtration (Whatman^® 42). The molybdate-reactive P (P_OI) was assessed using the Murphy and Riley colorimetric method adapted by Dick and Tabatabai (1977). This method has an accuracy of 95%–105% and a calibration curve with a range of 0.00 μg/ml to 2.50 μg/ml in an aliquot. To determine total carbonate-extractable (P_OT), a portion of the same extract was analyzed by ICP-OES (Varian Vista Pro CCD Simultaneous ICP-OES) using certified ICP-MS reference standards from SCP Scientific (Mississauga, Canada) for all desired elements. All phosphorus analysis performed was checked with a Standard Reference Material from the National Institute of Standards and Technology (NIST), United States. The standard used was the 2711a Montana II Soil. Based on the procedure used for ICP-OES analysis, the standard recovery was within 5%.

TABLE 1

TABLE 1. Summary of methods used to identify phosphorus pools and their symbols.

Total soil reactive P determined using acid persulfate digestion (P_SPT) provides results comparable to P_OT, but P is assessed by colorimetry rather than by ICP-OES (Kovar and Pierzynski, 2009). Analysis of P_SPT was performed using 0.5 g of potassium persulfate and 3 ml of 2.5 M H₂SO₄ with 5 ml of the P_OI extract. This mixture was boiled for 30 min, cooled and the pH adjusted with 2 M NaOH prior to assessing P by colorimetry, as for P_OI (Murphy and Riley, 1962; Dick and Tabatabai, 1977). This method was applied only to soil from the freeze-thaw experiment. The accuracy and calibrated standard curve procedure were the same as used above.

The concentration of total bioavailable elements [P_AD, Metal_AD] was determined using reverse aqua regia extraction (USEPA Method 3051A; Richardson and Reddy, 2013). One Gram of oven dried (105°C) and sieved soil (2 mm) was digested with 9 ml HNO₃ and 3 ml HCl in tightly sealed Teflon reaction vessels. Cold digestion for 24 h was followed by heated digestion at 100°C for 12 h. After extraction, samples were filtered (Whatman^® 42), diluted, and analyzed by ICP-OES (described above).

We used the chloroform fumigation extraction method to assess microbial biomass P (P_MBIO). This method measures P stored in microbial biomass by killing microbes with chloroform and measuring differences in P_OI extracted from fumigated and unfumigated soils (Carter and Gregorich, 2006). Soils from all sample sites were thoroughly mixed to make a composite sample, which was corrected for moisture content. Soil samples were placed in a desiccator and fumigated with boiled chloroform for 24 h. After incubation, the soil was extracted for P_OI, as described above. The soil P_MBIO was determined using the equation used and developed by Brookes et al. (1982), as modified by Voroney, Brookes, and Beyaert (Carter and Gregorich, 2006). By comparing unamended soil with soil spiked with 122.34 µg of P, the fraction of P_MBIO can be calculated using a Kep factor of 0.4, based on the assumption that 40% of P_MBIO is extractable by the chloroform method (Brookes et al., 1982; Carter and Gregorich, 2006). It is important to note that soils for the Olsen extractions (P_OI and P_OT) described above were air-dried, which can be expected to release P from microbial biomass. However, P_OI extracts captured on average just 35% of microbial biomass-P. Thus, we look at P_OI/P_OT and P_MBIO as distinct but related pools.

Soil microbial DNA extraction and sequencing

To identify seasonal shifts in microbial community structure, soil cores (2 cm diameter) were collected into sterile Falcon tubes during regular soil sampling in the wetland. One core was collected at each site for a total of nine cores per sampling day (Supplementary Table S1). Cores were transported to Guelph and immediately placed into a −80°C freezer until extraction for DNA, which was performed using the DNEasy Powersoil Pro Kit. After extraction, 16 S rRNA was quantified and stored at −80°C prior to submission to Metagenom Bio (Waterloo, ON, Canada) for library prep and 16 S rRNA sequencing on an Illumina MiSeq Machine.

Demultiplexed sequences were processed using the DADA2 pipeline (v1.8) (Callahan et al., 2016) in QIIME 2 (Caporaso et al., 2010; Bolyen et al., 2019). Reads were truncated at decreasing quality and quality filtered using default parameters (Figure S3, S4). Taxonomy was assigned using a naïve Bayesian classifier in QIIME2 with scikit learn. Taxonomy was assigned using SILVA release 138. Amplicon sequence variant (ASV) tables were imported into R for downstream analysis (R Core Team, 2021). All analyses were performed using the packages “phyloseq”, “vegan”, “zcompositions”, “microbiome”, and “ggplot2” (McMurdie and Holmes, 2013; Palarea-Albaladejo and Martín-Fernández, 2015; Wickham, 2016; Lahti and Shetty, 2017; Oksanen et al., 2018; Leite and Kuramae, 2020). Additional information on DNA extraction and sequencing conditions is provided in Supplementary Information.

Changes in soil P and microbial respiration during freeze-thaw cycling

Mesocosm incubation experiments

We investigated whether FTC changes the mobility and reactivity of soil P using acrylic boxes with stacked ports for water flow (Supplementary Figure S2). Each box had a volume of 2420 cm³. Soil collected from the upper 10 cm of wetland sites 1-9 was combined (2 kg each) and homogenized by hand. The soil (2280 g) was placed in the boxes at a field representative bulk density of ∼0.95 g/cm³. Boxes were insulated on sides and bottom with a 5 cm thickness of Styrofoam insulation to ensure that freezing occurred from the exposed upper surface down. There were three treatments: control, no-flow, and flow. Control mesocosms were held at 4°C throughout the duration of FTC cycling for the treated mesocosms and did not receive water flow. Deionized water was introduced using a peristaltic pump at a rate of 0.8 ml/min to slowly saturate the soil over 24 h, with an approximate pore volume of 1200 cm³. There were two replicates for each treatment, and each treatment was repeated three times, for a total of six replicates. The treatments consisted of incubation at −10°C for 10 days, followed by up to six, 48-h FTC at ± 4°C. Each experiment lasted 34 days. The no-flow mesocosms received only the initial saturating flow of water, followed by up to six freeze-thaw cycles. To assess soil P throughout the experiment, a 250 ml volume of soil was removed from one of six sections after 1, three and six FTC. No section was sampled more than once. The flow mesocosms received one pore volume of water (1200 cm³) via peristaltic pump after each 24-h thaw period. The discharged water was collected continuously in 250 ml fractions and assessed for DRP. Soil in the flow mesocosms was assessed only after six FTC to avoid changes to soil structure that would affect water flow. P_OI, P_SPT and P_MBIO were determined after 0, 1, 3, and six FTC treatments for the no flow treatments, and after six FTC for the flow treatments to avoid changes to soil structure that affect water flow. Dissolved P was analyzed in the discharged water from the flow treatment after each of the six thaw periods, after filtration (cellulose acetate 0.45 µm. Each of the three leachate samples collected per thaw were analyzed in four 10 ml replicates.

Control mesocosms comprised of soil that was not frozen were used to track the same fractions of P, designated as 0, 1-FTC-C, 3-FTC-C, and 6-FTC-C. All soil P fractions were quantified as described above for the field study.

Microbial respiration and P fractions after FTC

The impact of freeze-thaw cycling on soil P and microbial respiration was investigated using composited wetland soil (210 g) placed in small (250 ml) glass jars with screw-top lids and rubber seals. Soil from three samples dates (December 2020 April 2021, June 2021) was used in separate treatments and each treatment was done in triplicate. The jars were covered with Kimwipes secured by elastic bands to allow air exchange and were placed into Styrofoam blocks in an environmental chamber for temperature cycling. This facilitated the freezing of the soil from the exposed surface down. The soil in the jars was sampled destructively. Soil Olsen P was determined as described above. Soil respiration assays were conducted using the same soil in identical triplicate treatments. The soil was sampled after 0, 1, 3, and six FTC once the soil was thawed. Four portions of 15 g of soil from one jar were each placed into 250 ml glass serum bottles, sealed and crimped shut. These four vials represented time 0, 4-, 6-, and 24-h post-thaw, in order to track the recovery of respiration. The headspace of each vial was destructively sampled by drawing and analyzing three consecutive 5 ml gas samples with a syringe. The CO₂ concentrations were analyzed by gas chromatography, with a flame induction detector (FID) coupled to a methanizer (SRI 8610 GC). The standard curve used had detection points at 0, 5000, and 50,000 ppm. The retention curve for each sample (done in triplicate) was converted using the curve.

Statistical analyses

All statistical analyses were performed in R (R Core Team, 2021). To compare soil parameters over time (sampling date for field study and time point for lab studies), linear models were fit in R using REML (Bates et al., 2015; Wickham et al., 2021). For field samples, effects of sampling site on data were accounted for using Satterthwaite’s method and pairwise comparisons were made among dates using estimated marginal means (Lenth et al., 2021). For freeze-thaw study samples, linear mixed models were used and means were adjusted using the Tukey-Kramer adjustment, accounting for block as a random effect. Pearson correlations were used to determine correlations among variables. Principal component analysis (PCA) was performed to visualize trends in soil chemistry (Kassambara and Mundt, 2020), and trends were verified using multivariate analysis of variance (MANOVA). Soil P_MBIO was analyzed statistically using ANOVA as there was no random effects factor because the soil was homogenized across sample sites.

Respiration in the microcosms was compared across treatments at each hour of sampling. A linear mixed model fit using REML assessed the significance of CO₂ response to treatment and time. The model accounted for block as a random effect factor. The type III analysis of variance table with Satterthwaite’s method. To compare treatments within hours used an estimated marginal means adjusted to Sidak (Lenth et al., 2021). In a pairwise fashion, each treatment was compared to the other within each hour.

Results

Changes in environmental conditions prior to and during the study period

Southern Ontario typically experiences cold winters with freezing temperatures and intermittent warming. Air temperatures remained below freezing for most of the period from December 2020 until early March 2021, with an extended period in February when temperatures at night were −10°C or less. Despite these conditions, the soil appeared to never reach freezing temperatures beneath the snow layer, which was around 10 cm just before the onset of thaw conditions. The wetland is partly flooded each spring when the river overtops its banks during snow melt, which typically occurs in mid-to late-March. On 12 March 2021, the two transects closest to the river were flooded to a maximum depth of around 1 m. The treed site was not flooded.

The concentrations of dissolved reactive P were generally low (<0.1 mg/L) but were higher (0.2–0.4 mg/L) during the flood (data not shown). At this time the flood waters had encroached and covered the creek sampling location. This change reflects the widespread mobilization of P in the landscape; however, such elevated concentrations were only observed for a short time (<1 week) and the river and creek quickly returned to pre-flood concentrations. After the spring thaw and flood period, soil and air temperatures were coupled. The soil water content generally varied between 0.4–0.5 m³/m³, with lower values during late spring to summer and wetter conditions in fall and winter. Soil conductivity was lower in winter and early spring (0.8–1 mS/cm), and higher in summer (1–1.4 mS/cm), in agreement with the general pattern for soil moisture (data not shown).

Soil phosphorus pools during the spring to summer transition

There was both spatial and temporal variability in soil P pools in the landscape (Figure 2). Concentrations of P_OI in the wetland soil were consistently higher than in the area with trees (t (76) = 4.376, p < 0.001). This corresponds to the higher variability in vegetation, water flow and topography in the wetland. The predominance of standing mature trees may cause rapid scavenging of the most reactive form of P at low concentrations, or provide abundant living root biomass for P sorption in the rhizosphere. In addition, fewer sites were sampled in the treed zone. In the wetland, P_OI increased from early spring until May, corresponding to increased soil warming and microbial growth. A subsequent decline in June corresponded to the period of most rapid plant growth (Figure 2A). In contrast to P_OI, P_OT is more closely correlated between the treed zone and the wetland during the spring and is of similar magnitude (Figure 2A,B). The soil reactive P fractions made up a relatively small proportion of total acid extractable P (P_AD), around 2–5%.

FIGURE 2

FIGURE 2. Soil P in the most reactive fractions across the site from 2020–2021. (A) Soil P_OI, (B) Soil P_OT, and (C) P_MBIO (wetland only). ^abc represent significant differences among dates using estimated marginal means from the linear mixed model at p < 0.05 among wetland samples only. There were no significant differences among P_MBIO samples.

The trend for microbial biomass P in the wetland shows increasing values from early spring into May, followed by a large decrease in June. This suggests that the availability of highly reactive P is closely tied to the incorporation of P into microbial biomass. The data from the flood date appears to contradict this trend; however, the data is likely skewed because sites near the river were inaccessible on this date due to the high water level. The high DRP concentration during the flood period supports that the river delivered P to the wetland site, an ephemeral occurrence in the spring (Figure 2C). P_MBIO in the treed area followed a trend similar to that in the wetland, but the values were generally lower and less variable, with a range from 1.1 to 14.2 mg/kg during the spring. It is notable that P_MBIO was relatively high in December in the wetland, when both reactive P fractions were at their lowest values. This indicates that reactive P was immobilized within the microbial biomass at the beginning of winter. P_MBIO was also high relative to P_OT in May. This suggests that reactive P may be immobilized by microbes as the growing season progresses. As samples were composited across the site, there were no statistically significant differences among dates. In general, methods to determine P_MBIO result in large standard error (Bilyera et al., 2018), and the natural heterogeneity of soil also contributes substantial standard error (Xu et al., 2013).

We evaluated the chemical and physical data by date and by location, with Transect 1 adjacent to the river (sites 1–3) and Transect 3 furthest from the river (sites 7–9). The concentrations of most of the acid-extractable elements, averaged for all sites and dates and evaluated using linear mixed models, did not appear to vary significantly during the year (F (62) = 20.537, p < 0.001). Total calcium is the exception; Ca_AD was significantly lower (17,500 ± 3 107) during the flood event compared to the average for all other dates (50,910 ± 19,989). A closer look at the spatial distribution of the elements reveals, however, that calcium was particularly enriched in the sites proximate to the river, which could not be accessed when the river overtopped its banks (Table 2). Because the Grand River flows through Paleozoic limestone and more recent glacial deposits, it is not surprising that Ca is high where the influence of the river is the strongest (Gao, C., 2011).

TABLE 2

TABLE 2. Average total element concentrations December 2020 to June 2021 in the wetland across transects.

Similarly, Fe, Al and P were always more enriched at the sites that could be sampled during the flood, which most likely skewed the average to higher values on that sampling date. It is likely that more consistently saturated conditions and exposure to river flows through multiple flow paths caused the chemical reduction of FeIII and the consequent release of Al and P that were associated with the Fe oxides (Darke and Walbridge, 2000; Forsmann and Kjaergaard, 2014; Khan et al., 2022; Zhang et al., 2022).

Principal component analysis (PCA; Figure 3) helped to identify the trends in the physical and chemical data. Two major trends were clear in the PCA, highlighting temporal and spatial differences in the wetland. Along dimension 1 (representing 42.2% of the variation in the dataset), acid-extractable Al, Fe, and P are associated more strongly with Transect 2 and especially Transect 3, while total reactive P (P_OT) and acid-extractable Ca (Ca_AD) are associated more strongly with Transect 1. Along dimension 2 (representing 26.6% of the variation in the dataset), water content, microbial biomass P (P_MBIO) and inorganic reactive P (P_OI) are closely related, and similar across Dec-May (this is further validated by Figure 2, which shows a large dip in both P_OI and P_MBIO in June). In contrast, also along dimension 2, soil temperature appears to be inversely associated with water content, P_OI, and P_MBIO, reflecting a seasonal shift toward dryer conditions in June. It is notable that P_MBIO and temperature are inversely related, which could reflect the role of microbes as the initial traps in early spring for P that is transferred to plants as the cold soil warms and microbial biomass is mineralized.

FIGURE 3

FIGURE 3. Principal component analysis of key chemistry and soil sensor data from the wetland during the 2020–2021 field sampling season. Shapes indicate sampling day while colors indicate transects. Transect 1 = sites 1, 2, and three closest to the Grand River, transect 2 = sites 4, 5, and six in the middle, and transect 3 = sites 7, 8, and nine furthest from the river.

Microbial community dynamics

Microbial communities in the wetland were broadly dominated by common soil phyla, such as Actinobacteriota (13 ± 4%), Acidobacteriota (16 ± 4%), and Proteobacteria (31 ± 5%) (Figure 4). At the order level, the most abundant taxa were the Vicinamibacterales (12 ± 3%) of the Acidobacteria, and the Rhizobiales (10 ± 3%) and Burkholderiales (10 ± 3%) of the Proteobacteria. Microbial communities were significantly different between the two sequencing runs, possibly indicating differences in extraction and sequencing protocols (PERMANOVA: F (1,107) = 12.546, p = 0.001) (Supplementary Figures S3, S4). However, as noted above, soil moisture increased around October 2020, with a corresponding decrease in soil electrical conductivity; it is possible that the shift in microbial community structure reflects significantly different conditions (Supplementary Figure S5).

FIGURE 4

FIGURE 4. Microbial relative abundance and chemical parameters associated with microbial communities in the wetland samples. (A) Mean relative abundance of prokaryotic phyla representing >0.5% relative abundance in at least one sample. (B) Mean relative abundance of most prevalent prokaryotic orders across the 2021 growing season. (C) Redundancy analysis of clr-transformed prokaryotic ASVs (amplicon sequence variants) and selected environmental variables from December 2020-June 2021 in the wetland. AD. Fe/Al/Ca/PT are total concentrations from acid digestion.

The major notable difference between soil collected in October 2020, and in December 2020 and subsequent dates, was an increased relative abundance of the orders Nitrosophaerales (1 ± 1% in October and 12 ± 5% in December) and Rhizobiales (6 ± 4% in October and 14 ± 3% in December); Nitrosphaerales are ammonia-oxidizing archaea, and Rhizobiales are important bacteria in symbiotic N-fixation (Supplementary Figure S5). It is notable that these groups shifted with P_MBIO, which was highly seasonally dependent. During the growing season, the relative abundance of members of the order Rhizobiales tended to be higher in early and late spring (15 ± 3% in March 2021 and 14 ± 2% in June, with 12 ± 1% and 12 ± 3% in April and May, respectively) (Figure 4B); this follows a pattern nearly opposite to P_MBIO. In contrast, the relative abundance of members of the Nitrososphaerales followed a pattern similar to P_MBIO, with a trend to lower values in March and June (10 ± 4% and 12 ± 5%) compared to April and May (16 ± 3% and 14 ± 3%).

Downstream analyses were performed using data from December 2020-June 2021 only, as those have more complete geochemical metadata. Microbial community structure, like P_MBIO, was significantly associated with date of sampling (PERMANOVA: F (5,36) = 1.3151, p = 0.004). PERMANOVA analysis also revealed significant differences in microbial community structure among areas of the site (F (8,33) = 1.4376, p = 0.001). The analysis showed that P_AD (F (1,40) = 2.5568, p = 0.001), Fe_AD (F (1,40) = 2.4771, p = 0.001), Al_AD (F (1,40) = 2.4797, p = 0.001), and Ca_AD (F (1,40) = 2.8911, p = 0.001) were the only soil chemistry parameters associated with differences in microbial community structure over time. These parameters were, therefore, selected for redundancy analysis (Figure 4) to correlate shifts in microbial communities with environmental variables (Figure 4). 44.54% of variance in microbial communities was explained by RDA axis 1 and 22.59% by axis 2. Similar to the PCA (Figure 3) the RDA shows a negative relationship between Ca_AD and P_AD/Al_AD/Fe_AD (Figure 4). This suggests that soil from Transect 1, which had relatively high total Ca and low total P, Al and Fe, had lower microbial diversity than that from Transects 2 and 3, which experienced less saturated conditions and less direct contact with the river.

Changes in soil phosphorus and respiration during freeze-thaw cycling

As observed (Figure 5), despite extended sub-zero temperatures, the soil did not freeze. To explore one future climate scenario in which soil temperature fluctuates in early spring increase due to reduced snowpack (UNFCCC, 2008; Kieta et al., 2018), we conducted studies using mesocosms to assess significant changes in soil P fractions in response to freeze-thaw cycling. There were two treatments: 1) FTC only; and 2) FTC with low water flow between periods of freezing to mimic seepage of pore water. The control consisted of the same soil but kept at 4°C throughout. In the treatments, three or fewer FTC did not affect P_SPT, while P_SPT was significantly higher after six FTC in the treatment without flow compared to the control, but not compared to treatments with fewer FTC (Table 3; Figures 6, 7). Fluctuations in P_OI followed the same trend but the differences were lower in magnitude. There was little change in P_MBIO with FTC treatment, although P_MBIO increased over time in the control treatments. In the flow mesocosm treatments, water was introduced at a rate to mimic seepage through the soil between freeze-thaw cycles. The leachate collected after the first thaw in the flow mesocosms had the highest concentration of DRP, as expected due to displacement of pore water P (Figure 7). Each successive FTC resulted in lower or stable values for DRP concentration, until six FTC had taken place. It is notable that DRP increased after 6FTC despite the dilution of pore water P by the inflow water after successive FTC, and in agreement with a significant increase in P_OI and P_SPT. This reflects the marked increase that we observed for P_OI after six FTC for both flow and non-flow treatments and suggests that reactive soil P is mobilized by multiple cycles of freezing and thawing despite the potential for increased P adsorption due to the impacts of freezing on soil structure, as observed by others (Wang et al., 2007).

FIGURE 5

FIGURE 5. Air and soil temperature data (°C) from the ATMOS-14 and TEROS 11/12 sensors on the Zentra data logger between December 2020 and June 2021. Circles indicate sampling dates over the field season. Temperatures were monitored from June 2020–December 2021 but are not shown in the figure.

TABLE 3

TABLE 3. Microcosm P _OT and P _OI after 1, 3, or six FTC and in control soil kept at 4°C after the concurrent FTC treatments (FTC-C).

FIGURE 6

FIGURE 6. Soil P pools after freeze thaw cycling in mesocosms. (A) Soil P_OI in mg/kg (B) Soil P_SPT in mg/kg. (C) Soil P_MBIO in mg/kg. “0” refers to samples taken prior to chamber exposure, and “FTC” indicates how many freeze-thaw cycles were performed prior to sampling. “FTC-C” refers to control soils kept at 4°C after the concurrent FTC treatments. ^abc represent significant differences among dates using estimated marginal means from a linear fixed-effects model at p < 0.05. There were no significant differences in P_MBIO.

FIGURE 7

FIGURE 7. Leachate DRP (mg/L) following FTC treatment. ^abc represent significant differences among dates using estimated marginal means from a linear fixed-effects model at p < 0.05.

Changes in soil P and soil respiration during freeze-thaw cycling

To relate P dynamics to microbial activity in response to FTC, we conducted additional freeze-thaw experiments in 250-ml jars that allowed us to evaluate respiration and shifts in the P fractions together. As observed for the mesocosm studies, soil P_OI varied significantly only after six FTC (Table 3). Values for P_OI were 373% higher after six FTC, compared to the untreated soil, and 290% higher for P_SPT, which mirrors the trends observed for the mesocosm treatment. Although the increase in P_SPT was highest after six FTC, there was a significant increase after three FTC; this likely reflects a higher impact of freezing processes in the soil, similar to the observation in the mesocosm study. The control treatments generally had lower or comparable concentrations of reactive soil P. The exception is P_SPT after three cycles (3 FTC-C), which was higher and similar to six FTC-C, indicating that P_SPT was stabilized (Table 3). Freeze-thaw cycling stimulated net CO₂ production compared to the control and to the untreated soil at all treatment levels (Figure 8). Respiration was suppressed after six FTC treatments compared to one and three FTC.

FIGURE 8

FIGURE 8. Clockwise (A) Soil respiration (ppm/g) over the first 24 h pre-incubation; (B) after 1FTC (C) after three FTC; (D) after six FTC.

Discussion

The highest concentrations of dissolved phosphorus in stream water were observed during spring thaw, when the river overtopped its banks. This coincided with peak flow conditions. Although P fluxes in streamflow were not measured in this study, such periods are typically when the majority of annual P losses are observed in a given year in the study region (e.g., Macrae et al., 2007; Macrae et al., 2010; Irvine et al., 2019) as well as other temperate landscapes (e.g., King et al., 2015; Plach et al., 2019). Soil inorganic reactive P (P_OI) was also high at this time, while total reactive P and microbial biomass P were low. This suggests that mineralization of P during the winter months increased P_OI at the expense of organic and microbial biomass P (P_MBIO). The timing of spring thaw and runoff is, therefore, critical for loss of reactive P to the downstream environment.

The seasonal changes in total reactive P follow a similar pattern in the wetland and wooded sites. The response of microbial biomass P is, however, more dynamic in the wetland and is closely related to P_OI. The soil P_MBIO fraction is well understood in terms of function and size (Richardson and Simpson, 2011). The concentration of P_MBIO in the upper 30 cm of wetland soil globally has been estimated at around 71 mg kg⁻¹ with a range of 0.6–345 mg kg⁻¹ (Wang et al., 2022), in agreement with P_MBIO determined for the soil at our site. How this fraction responds within riparian systems to environmental change remains poorly understood. The P_MBIO regulates P_OI availability in soil environments especially during the non-growing season when there is little to no plant uptake (Richardson and Simpson, 2011; Sorensen et al., 2018), as we inferred for the period of seasonal thaw in the wetland. In the wooded area, with the dominance of mature deciduous trees and fewer seasonal understory plants, microbial biomass P is relatively stable. In the wetland, the increase in microbial biomass P in spring after the freshet and before the emergence of plants was followed by a rapid decline that corresponded with the establishment of dense vegetation. The tendency for total reactive P to increase during the growing season suggests that some microbial P is transferred into this pool as spring progresses. Although P_MBIO is a relatively reactive form of P, it is distinct from total Olsen-P. The return to higher concentrations of P_MBIO and the decrease in P_OT by December indicates the return of phosphorus to the microbial biomass pool during seasonal plant senescence. The relationship between P_MBIO, water content and orthophosphate (P_OI) supports that the mobilization and uptake of P corresponds to conditions that foster microbial proliferation. It is notable that while P_OI and P_MBIO appeared to be seasonally associated with soil temperature and water content, P_OT and P_AD were more spatially dependent, with higher P_OT and lower P_AD in Transect 1 than in Transects 2 and 3. The soil conditions may be more similar across the wetland when river levels decline in late spring, although it is difficult to know with certainty that P_MBIO is spatially dependent because samples were composited across the site for each given date. The chemical forms of soil organic P are highly diverse and are poorly understood compared to inorganic P (Huang et al., 2017), which prevents a more detailed interpretation of P cycling at our site.

Rivers in southern Ontario are high in calcium, which can explain the relatively high concentrations of Ca in the wetland at the river’s edge. The soil pH in the wetland averaged 7.5, with little variation across sites. This suggests that Ca should control P solubility (Lindsay et al., 1989; Dunne and Reddy, 2005). We observed, however, a close correlation between Ca and P only for the transect adjacent to the river, and only for reactive P. The iron mottles in the soil that occurred mainly in Transects 2 and 3 support that Fe redox cycling is seasonally active. This indicates that there may be a metal oxide trap for P despite bulk soil conditions that favor precipitation with calcium rather than reaction with Fe and Al oxides (Lindsay et al., 1989). We speculate that phosphorus pools are likely controlled at least in part by microsites that are chemically distinct from the bulk soil through redox dynamics regulated by microbial activity. Such sites can occur in soil concretions and nodules (Gasparatos et al., 2019) and on the roots of wetland plants, where Fe oxide plaques are known to trap P (Xu et al., 2009; Zhu et al., 2018). More research is, however, needed to confirm whether metal oxides play a role in seasonal sequestration of P in this wetland. There is less potential for P to be immobilized by mineral associations in anoxic conditions; although the reduction of FeIII in iron oxides to FeII in the presence of phosphate can result in vivianite, this mineral is highly unstable in oxidizing conditions (Lindsay et al., 1989; Glasauer et al., 2003). Microbial community structure was broadly correlated with Ca_AD, Al_AD, Fe_AD, and P_AD, indicating broad shifts in taxa in response to elements that are related to redox processes. There is likely also a spatial component, with microbial communities responding to harsher conditions on the front edge of the wetland, in particular less oxygen and higher leaching of nutrients. Interestingly, despite the strong correlation between acid digestible P and Fe, major taxa associated with Fe oxidation and reduction (ex. Thiobacillus, Leptospirillum, Geobacter) represented low overall relative abundance (<1%) in our study.

We also observed that microbial communities dominated by N-fixing bacteria (members of the order Rhizobiales) shifted to those dominated by ammonia oxidizing archaea (the order Nitrososphaerales) over the course of the growing season. In general, soil N-mineralization and nitrification rates tend to increase as the soil warms. In the early spring, the oxidation of ammonia during the seasonal increase in mineralization is likely, similar to the correlation between P_OI and P_MBIO that we observed during the spring (Zhang et al., 2008; Contosta et al., 2011). Other studies have found correlations between soil P and abundances of ammonia oxidizing archaea (Tang et al., 2016), signalling shifts toward increased N mineralization as microbes utilize available inorganic soil P. Our study is unique because we included microbial biomass P, providing a more direct insight into the ability of microbes to assimilate P over the growing season.

Agricultural areas are responsible for large losses of P in runoff, and wetlands can help to attenuate such losses to downstream waters (Osborne and Kovacic, 1993; Dunne and Reddy, 2005). The spatial data from our study suggests that the Grand River wetland is more effective at sequestering P with distance from the river (Cooper and Gilliam, 1987); however, the soil situated further from the river likely intercepts higher amounts of P from the upland soil and vegetation. The significantly lower concentrations of acid extractable Fe, Al and P adjacent to the river suggest that leaching is enhanced by more saturated conditions and by contact with flowing water. The high concentrations of DRP in the creek and river during annual spring flooding indicates widespread losses of P broadly across the agricultural landscape over a relatively short period of time. Although floodplains can serve as traps for P during seasonal flooding (Gillespie et al., 2018), the gully riparian wetlands that are typical for the hilly environment of the upper Grand River watershed are small relative to the upland areas that generate runoff. One potential solution to reduce P inflows to the river is to divert surface runoff that currently flows to the creek to areas designed to attenuate dissolved and suspended P loads (Dunne and Reddy, 2005).

We were interested in whether freeze-thaw cycling in the wetland would impact reactive soil P, in order to predict whether anticipated changes in soil freezing for cold temperate soils are likely to increase the losses of P to surface waters. According to the buried temperature sensors, the wetland soil did not freeze despite an extended cold period. This suggests that FTC was not the driver of changes in P fractions during the year of our study. It is possible that freezing and the mobilization of dissolved reactive P occurred above the 10 cm depth of the soil probes. However, the probes were placed at the depth of soil sampling. The soil surface beneath the thin snow and ice cover was always wet during site visits. It has been suggested that where FTC occurs, it is most likely that there will be multiple cycles over the course of the season, such as in early spring (Henry, 2007). In our laboratory study, we observed that multiple (>3) FTC were required to increase reactive soil P. A decrease in dissolved P with FTC cycling may reflect increased adsorption of P to disrupted soil aggregates (Wang et al., 2007; Yu et al., 2011; Li et al., 2022) and a consequent increase in reactive soil P; freeze-thaw cycles have cumulative impacts on the disruption of soil aggregates (Fitzhugh et al., 2001). It is difficult to draw broad conclusions on how FTC impacts the size of P pools in wetland soil due to the small number of studies and the substantial experimental differences between them, which persist in part because of the difficulties in isolating a single process that can explain the impact of soil freezing on phosphorus (Vaz et al., 1994) Our mesocosm studies were designed to mimic natural conditions as closely as possible (Henry, 2007), given experimental limitations. Nevertheless, microbial biomass is more readily disrupted by freezing if microbes are not adapted to freezing conditions, which was the case for the wetland soil used in our study (Grogan et al., 2004; Bölter et al., 2005).

Our results indicate that FTC decreases microbial biomass P in wetland soil only after multiple cycles. Because FTC increased the reactive P fractions, we infer that P_MBIO was transferred into these fractions. This represents the first report of this phenomenon. More broadly, the loss of total microbial biomass due to increased freeze-thaw cycling has been suggested to offset the impacts of soil warming by decreasing the activity of soil microbes even during expected periods of high activity subsequent to FTC (Sorensen et al., 2018). This interesting result suggests that an increased return of available P to soil from biomass during FTC is not necessarily sufficient to offset the negative impacts of FTC. Others have also observed that microbial respiration is stimulated by FTC (Schimel and Clein, 1996). The burst of CO₂ that we noted upon thawing may be due to the release of CO₂ that was trapped in the soil during freezing, among other experimental artifacts (Song et al., 2017). Nevertheless, the marked difference in respiration after six FTC corresponds to the significant changes in the P fractions after six FTC, and agrees with other research showing disruption of microbial cells and decreased respiration with FTC for a wide range of soils (Song et al., 2017). Our results support that soil freeze-thaw cycles have the potential to modify soil P dynamics.

The importance of wetlands in riparian zones within agricultural watersheds for retaining dissolved and particulate nutrients is well known under warm conditions (Gillespie et al., 2018), but less so during cold conditions (Kieta et al., 2018). The necessary wetland physical parameters for retaining soil P, such as size, slope and elevation with respect to the adjacent water body, are poorly defined under any conditions. There is also limited information on groundwater-surface water interactions in and across riparian areas, and if and how seasonal shifts in flow paths may relate to the mobilization of P from these zones. Understanding the soil processes that control seasonal losses of phosphorus from riparian wetlands can help to predict conditions which may overwhelm the natural capacity of the soil to sequester P.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, BioProject Accession: PRJNA849993.

Author contributions

JC and SG designed the study, and SG supervised. MM provided guidance throughout the project. JC performed all sample collection and lab work, while KM and JC collaborated on data analysis. KM performed all microbial data analyses and interpretation. JC wrote the first draft of the manuscript and all authors contributed to subsequent drafts. This manuscript has been approved by all authors.

Funding

This research was funded by the Ontario Ministry of Agriculture, Food and Rural Affairs (OMAFRA) Grant #UofG 2018-3254.

Acknowledgments

We thank Dr. Al Mattes of NatureWorks for providing site access and support for the field work, and Dr. James Longstaff, Dr. Paul Voroney, Peter Smith and Joseph Radford at the University of Guelph for advice and technical assistance.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fenvs.2022.983129/full#supplementary-material

References

Bates, D., Mächler, M., Bolker, B. M., and Walker, S. C. (2015). Fitting linear mixed-effects models using lme4. J. Stat. Soft. 67. doi:10.18637/jss.v067.i01