Phosphate oxygen isotope ratios in vegetated riparian buffer strip soils

The oxygen isotopic composition of phosphate (δ18OP) in soils and surface water bodies has been used to trace terrestrial P inputs into aquatic ecosystems. However, enhanced biological activity in vegetated riparian buffer strips (VBSs) may lead to an alteration of δ18OP values. The objective of this study was to assess whether enhanced biological cycling of P in VBS soils can be identified using δ18OP values. For this purpose, we sampled temperate grassland soils at various depths along a VBS to grassland transect. Here, we combined sequential P soil extracts with an analysis of δ18OP values. Soil P pool concentrations tended to decrease significantly along the transect from the VBS to the grassland soils; the strength of this relationship varied with P extract, sample depth, and inorganic or organic bonding form. For the δ18OP values of the 1 M HCl‐extractable P we observed a significant negative trend along the VBS to grassland transect, indicating a tendency for accelerated rates of biological cycling of P within the VBS soil profile compared with the upslope soils. We conclude that oxygen isotope‐based assessments of P source contributions to freshwater bodies should consider the enhanced biological turnover of P in VBS soils.


Vadose Zone Journal
A reduction of nutrient inputs to surface waters can be achieved by establishing vegetated strips on riverbanks, called vegetated riparian buffer strips (VBS) (Lee et al., 2003;W. Roberts et al., 2012;Vidon et al., 2010). Hence, VBS have been identified as ecosystem control points for nutrient retention and cycling (Bernhardt et al., 2017;Hedin et al., 1998;McClain et al., 2003). Vegetated riparian buffer strips physically filter soil particles from water in overland flow, and thereby, reduce nutrient loads and their transport rates into surface waters, although the dynamics of transport and release processes of P from adjacent agricultural lands through VBS into streams is still poorly understood (Young et al, 2021). With the reduction of transport rates, that is, longer residence time of water, presumably also the infiltration of water and nutrients into VBS soils increases compared with nonvegetated riverbanks. Therefore, P transported with overland flow is not directly transferred to surface waters but is retained in soil geochemical (e.g., via sorption and precipitation processes) and biological pools (e.g., via uptake by roots and microorganisms) (Dosskey et al., 2010;S. Gu et al., 2017;Saia et al., 2021;Zhang et al., 2010).
The water table in riparian zones is often shallow, which leads to the submersion of deeper soil layers in water over prolonged periods of time. Given a sufficient supply of soil organic matter and electron acceptors, these submerged soil layers can develop a very low redox potential, which eventually leads to the reductive dissolution of ferric iron (Fe 3+ ) and the simultaneous mobilization of iron (Fe)-bound phosphate (Bonneville et al., 2004;Parson et al., 2017;Schwertmann, 1991). At the same time, biological activity may be elevated in soils of VBS, which can accelerate the rate of transformation of organic matter and nutrients (Jacinthe et al., 2015;Stutter et al., 2009;Tufekcioglu et al., 2001) and thus also of biological P cycling.
Biological P transformation in soils can be identified using the O isotopic composition of phosphate (δ 18 O P ) (Middelboe & Saaby, 1998;Tamburini et al., 2012). In the absence of biological activity, the δ 18 Ο P values do not change at temperatures below 70˚C (O'Neil et al., 2003). At ambient temperatures, only enzymatic activity alters δ 18 Ο P values by catalyzing the exchange of O between phosphate and water. The hydrolysis of organic P compounds by extracellular phosphatases leads to the incorporation of one (phosphomonoesterases) or two (phosphodiesterases) O atoms from water into phosphate. This O atom exchange involves a substratedependent kinetic isotope fractionation, which has been quantified for alkaline phosphatases (Liang & Blake, 2006), acid phosphatases (von Sperber et al., 2014), phytases (von Sperber et al., 2015;Wu et al., 2015), and phosphodiesterases (Liang & Blake, 2009). When plants and microorganisms take up a phosphate molecule from the soil solution, the ubiquitous intracellular enzyme pyrophosphatase will induce repeated exchange of O atoms between phosphate and water. Over

Core Ideas
• Vegetated riparian buffer strips accelerate biological P cycling. • Phosphate oxygen isotope values decreased along a vegetated riparian buffer strips to grassland transect. • Biogeochemical P cycling varies with soil depth and distance from a stream. • Cycling-related changes in isotope values complicate stream P source identification.
time, sustained biological cycling of P within organismic cells results in an equilibration of O isotopes between phosphate and water, which involves a temperature-dependent isotope fractionation between O in water and O in phosphate (Chang & Blake, 2015). Incubation experiments with algae and different strains of bacteria have shown that the isotopic equilibration occurs within a time span of hours to days (Blake et al., 2005;Paytan et al., 2002;Stout et al., 2014). Several studies have compared δ 18 Ο P values of sequentially extracted soil P fractions to differentiate abiotic and biotic processes in natural systems such as those along climatic gradients (Angert et al., 2011;Helfenstein et al., 2018), along soil chronosequences (K. Roberts et al., 2015;Tamburini et al., 2012) or within soil profiles at individual sites (Amelung et al., 2015;Bauke et al., 2018;Joshi et al., 2016). Further, δ 18 Ο P values have been used to identify sources of P input to freshwater bodies by analyzing δ 18 Ο P in water and in potential sources such as fertilizer, waste water, and upstream soils within the respective catchment (Elsbury et al., 2009;Ishida et al., 2019;Paytan et al., 2017;Pistocchi et al., 2017;Tonderski et al., 2017). However, no study has considered how δ 18 Ο P values are affected by the biogeochemical cycling of P in VBS soils, a crucial ecosystem control point on the interface of the terrestrial and aquatic P cycles. Therefore, the objective of this pilot study was to obtain first insights into P cycling in riparian soils, by assessing soil P concentrations and δ 18 Ο P values along a first-order stream VBS to grassland transect. As riparian ecosystems are subject to frequent fluctuations in soil moisture status, oxygen isotope analyses were carried out only for the HCl-extractable pool, which may incorporate a longerterm signature of biological cycling (Angert et al., 2012;Joshi et al., 2016;Helfenstein et al., 2018) as opposed to more soluble or microbial pools, which are highly susceptible to specific moisture and redox conditions on the day of sampling. We sampled soil profiles to 100-cm depth to include a gradient in redox conditions, from strongly reducing conditions close to the riverbank and in the deeper subsoils to oxidizing conditions further away from the riverbank and in the top-soils. We hypothesized that (a) due to the retarding effect of the VBS on P transport by overland flow, P will be accumulated at locations close to the stream across all P pools and (b) as a consequence, also the biological cycling will be enhanced as indicated by δ 18 Ο P values closer to equilibrium.

Study site and sampling
We conducted our study along the first-order stream "Eichenbach" at the University of Bonn Research Farm "Frankenforst" (7˚12′ 22″ E, 50˚42′ 49″ N), Germany, where soils developed from loess loam of > > 1-m depth. Mean annual precipitation is ∼700 mm and mean annual temperature is 9.0˚C. The VBS has a width of ∼3 m and is vegetated mostly by willow trees (Salix alba L.) and diverse bushes including hazelnut (Corylus avellana L.), blackberry (Rubus plicatus L.), and stinging nettles (Urtica dioica L.). The sampled area adjacent to the VBS is permanent grassland with 2-5% slope, which has been used as a cattle pasture for more than 25 yr. The VBS is protected from cattle access by a fence. A map and pictures of the field site, the sampling transect and a soil profile are given in Supplemental Figure S1. We sampled soils with a soil auger along a transect with sampling points located 1, 17, 20, 43, 49 and 55 m away from the stream following the local topographic slope. Each sampling point was sampled once and the obtained material was separated into four depths (0-10, 10-30, 30-50, and 80-100 cm, see Supplemental Table S1). The soil material from each depth and sampling point was then homogenized. Soils along the transect represent a transition from Fluvents close to the stream to Orthents for the most remote grassland sites (Gleysols to Cambisols according to WRB [2015]). We selected sampling locations with respect to soil horizons and pedogenic features indicating groundwater influence and redox potentials. The sampling locations and depths represent different oxic and anoxic soil horizons with the goal to encompass the whole soil moisture and redox gradient, from strongly reducing conditions close to the riverbank and in the deeper subsoils (Gleysol close to the riverbank) to oxidizing conditions further away from the riverbank (Cambisol at the most distant site). We determined horizons visually from soil cores and with the information from an adjacent soil pit. We collected samples and stored them field fresh at 4˚C. For each sample, we determined gravimetric water content of the field moist samples by oven-drying of a sample aliquot at 105˚C. We also determined total C, N, and S contents by dry combustion of a sample aliquot (vario MICRO cube, Elementar, Hanau, Germany). We analyzed gravimetric water content, C, N, and S soil content in duplicates. We did not record redox potentials at the time of sampling, because these are known to vary throughout the season (Mansfeldt, 2003). Instead, we focused on assessing an integrated signal of soil P cycling across several weeks to years. Based on the soil profile classification, the Cr layers can be considered as permanently anoxic horizons, while Bl and Bg horizons are only seasonally anoxic during periods of high water tables. All other soil horizons are considered to be oxic (Supplemental Figure S1).

Sequential P extraction
We extracted soils at 1:60 soil/solution ratio (w/v) following an adjusted sequential P extraction scheme to quantify operationally defined soil P pools (Tiessen & Moir, 1993).
In brief, we first extracted 0.5 g of soil with H 2 O to extract directly bioavailable inorganic phosphorus (P i ) from the soil solution. In the following step, the soil was shaken in a 0.5 M NaHCO 3 solution to assess moderately available P i and easily mineralizable organic phosphorus (P o ). Next, we extracted samples with 0.1 M NaOH, which targets less available P i fractions bound to Fe or Al and less available P o fractions, which are still involved in short-term transformations of P in soils (Cross & Schlesinger, 1995). After the alkaline reagents, we extracted samples using 1 M HCl, which is assumed to extract mostly stable Ca-bound P, for example apatite (Cross & Schlesinger, 1995). For each P extract, we determined P i concentrations on a spectrophotometer at 900 nm wavelength (SPECORD 205, Analytik Jena) using the molybdenum (Mo)blue method (Murphy & Riley, 1962). All extractions for P i analyses were performed in triplicate for each transect distance and depth (n = 3 per depth and transect distance, with a total of n = 72 per extract). Total P analysis was conducted only once for each transect distance and depth (n = 1 per depth and transect distance, with a total of n = 24 per extract) by inductively coupled plasma optical emission spectroscopy (ICP-OES; HORIBA 2 Jobin Yvon). We calculated P o as the difference between total P and the mean of the P i triplicates for each extract, distance, and depth. This resulted in n = 24 per P o extract.

Oxygen isotope analysis in phosphate
In preparation of the analysis of the oxygen isotopic composition of phosphate (δ 18 O P ), we sequentially extracted 30 g of each sample with 300 ml of 0.5 M NaHCO 3 , 0.1 M NaOH, and 1 M HCl and with 16-h extraction time for each step. The extraction steps with alkaline solutions were mainly used to remove organic matter from the extracts, only the final 1 M HCl extract was processed further for oxygen isotope analysis (Amelung et al., 2015). For purification of the 1 M HCl extracts, we followed the procedure of Tamburini et al. (2010) with precipitation, washing and dissolution of phosphate as ammonium-phosphomolybdate and magnesium ammonium-phosphate, followed by rinsing of the solute on a cation exchange resin (Dowex 50 × 8, 200-400 mesh, Sigma) and final precipitation of phosphate as silver phosphate (Ag 3 PO 4 ). Henceforth, we use δ 18 O P when referring to the δ 18 O P value of the 1 M HCl extract. We then analyzed the Ag 3 PO 4 precipitate of each sample for δ 18 O P in duplicates (i.e., two repeated laboratory measurements per distance and depth, resulting in a total of n = 48) using a Vario PYRO Cube (Elementar) with a C-based reactor connected in continuous flow to an Isoprime 100 isotope ratio mass spectrometer (IRMS; Isoprime). We set the furnace temperature of the reactor to 1,450˚C and resulting CO gases were concentrated by a purge and trap chromatography system. We calibrated results against an internal Ag 3 PO 4 standard (Acros Organics; δ 18 O = 14.2‰) and two benzoic acid standards (International Atomic Energy Agency; IAEA 601: δ 18 O = 23.1‰ and IAEA 602: δ 18 O = 71.3‰). We reported results in the conventional delta notation relative to the Vienna Standard Mean Oceanic Water (VSMOW). In order to assess the purity of the produced samples, we compared their oxygen yields with those of the Ag 3 PO 4 standards (Supplemental Figure S2).
We calculated theoretical equilibrium values resulting from the complete equilibration of oxygen isotopes between phosphate and water mediated by pyrophosphatase enzymes based on the equation (Equation 1) of Chang and Blake (2015): where T is the ambient soil temperature in [K] and δ 18 O P and δ 18 O W are the standardized 18 O: 16 O ratios of phosphate and water, respectively (reported in delta notation relative to VSMOW). Because δ 18 O P values of the HCl-P fraction integrate the isotopic signal over long periods of time (Angert et al., 2011;Joshi et al., 2016) (Allison & Barnes, 1983;Hsieh et al., 1998). Corresponding long-term air temperatures measured at a weather station of the German Weather Service at Cologne-Bonn airport, approximately 18 km North Northwest from the field site, were 13.7˚C. As a result, the expected δ 18 O P value for biologically cycled phosphate was 18.3‰ (±2‰ to again account for the evaporation-based uncertainty in soil water isotope values).

Statistical analyses
We conducted all statistical analyses in R (version 3.6.2; R Core Team). We used the lm function, which is a builtin R function, to implement the P pool linear models and test whether P i and P o pools from the four extracts (i.e., H 2 O, NaHCO 3 , NaOH, and HCl) varied with transect distance and depth. We used the lmer R function in the lme4 R package (Bates et al., 2015) to implement the δ 18 O P linear mixed effects models and test whether δ 18 O P varied with transect distance and depth. Specifically, we used a linear mixed effects model because each δ 18 O P sample was analyzed two times on the IRMS (i.e., repeated laboratory measurements) so it was necessary to account for the variation between samples as well as the variation with repeated IRMS measurements. All models were run with the full number of data available, that is, n = 144 for P i , n = 24 for P o , and n = 48 for δ 18 O P . In our analyses, δ 18 O P values and concentrations of each P i and P o pool extract represented the response variables. We used both transect distance and depth as predictors and compared model performance for three linear model cases: (a) distance; (b) distance and depth; and (c) distance, depth, and the interaction of the distance and depth. See Supplemental Tables S2-S4 for full descriptions of each linear model. In some cases, we log-transformed the response variable to ensure response variables were normally distributed.
We compared the performance of the three regression model cases using the Akaike Information Criteria (AIC). A lower AIC value indicates a better model fit and AIC scores that are more than two units apart are significantly different. Best fit linear models are highlighted in italics in Supplemental Tables S2-S4. Additionally, for each best fit model, we checked that model residuals were normally distributed and homoscedastic (i.e., had constant variance). We tested model significance at alpha = .05 and used R 2 to measure the percentage variation explained by the predictor variables included in the model.
Because the sampled depth intervals may not be broadly relevant beyond this field site, we also tested whether distance and the long-term redox conditions, as defined by the soil horizons, could explain patterns in δ 18 O P values. Here, redox condition was a categorical value in the model and was defined by whether the sample originated from the Ah layer (i.e., oxidized) or the Cr layer (i.e., reduced) of transect distances 1 to 43 m (distances 49 and 55 m were not influenced as strongly by the groundwater and therefore did not have a Cr horizon; Supplemental Figure S3). Similar to the methods described previously, we compared the performance of three linear mixed effects model cases: (a) distance; (b) distance and redox condition; and (c) distance, redox condition, and the interaction of the distance and redox condition. As before, we used a linear mixed effects model because each δ 18 O P sample was analyzed two times on the IRMS (i.e., repeated laboratory Vadose Zone Journal F I G U R E 1 Soil inorganic phosphorus (P i ) concentrations for the H 2 O, NaHCO 3 , NaOH, and HCl extracts vs. transect distance (n = 72 per extraction). Black line represents the median, top of box represents 75th percentile, bottom of box represents 25th percentile and points indicate individual observations. For linear model results associated with this figure, see Supplemental Table S2 measurements) so it was necessary to account for the variation between samples as well as the variation with repeated IRMS measurements. See Supplemental Table S5 for full descriptions of each linear mixed effects model.

RESULTS
A summary of gravimetric soil water content on the day of sampling as well as soil C, N, and S contents are given in Supplemental Table S1. Considering the complete soil profile, concentrations of all P i extracts were significantly greater at locations close to the stream compared with those of upslope grassland soils (Figure 1; Supplemental Table S2). However, transect distance alone explained 12 and 35% of the variance in P i concentrations for H 2 O and NaHCO 3 extracts, respectively (Supplemental Table S2). For P i concentrations of NaHCO 3 , NaOH, and HCl extracts, including the interaction of transect distance and depth explained 57, 73, and 42% of variation in the data, respectively (Supplemental Table S2). Including the interaction term between distance and depth, significantly improved our ability to explain the variation in P i concentrations for these extracts compared with a model that only included the additive effect of distance and depth (e.g., Model no. 13 in Supplemental Table S2). Therefore, P i concentrations of NaHCO 3 , NaOH, and HCl extracts depend on and vary with both transect distance and soil depth. The P o concentrations showed a similar general pattern as compared to P i concentrations, with significantly larger P o concentrations near the stream than in the upslope grassland soil (Figure 2; Supplemental Table S3). However, again we observed that models of P o concentrations from all four Vadose Zone Journal F I G U R E 2 Soil organic phosphorus (P o ) concentrations for the H 2 O, NaHCO 3 , NaOH, and HCl extracts vs. transect distance (n = 24 per extraction). Black line represents the median, top of box represents 75th percentile, bottom of box represents 25th percentile and points indicate individual observations. For linear model results associated with this figure, see Supplemental Table S3 extracts vs. transect distance improved when we included sampling depth as additional factor or within an interaction term, explaining 38-53% of data variation (Supplemental Table S3).
The δ 18 O P values of 1 M HCl-extractable P in topsoil and subsoil close to the stream were within the equilibrium range across the complete profile until 100-cm depth. With increasing distance from the stream, however, the δ 18 O P values increasingly deviated from the equilibrium, first in the subsoil and then further away from the stream also in the topsoil. The δ 18 O P values in the soils most distant from the stream differed by almost 10‰ relative to the ones found closest to the stream. Indeed, we observed a significant negative trend of δ 18 O P values with increasing distance from the stream (Figure 3; Supplemental Table S4). Sampling depth and transect distance alone accounted for 13 and 64% of vari-ation in δ 18 O P values, respectively (Supplemental Table S4). However, when we included the additive effects of sampling depth and transect distance in the linear mixed effects model, we captured 76% of variation in δ 18 O P (Model no. 3, Supplemental Table S4).

Phosphorus retention in VBS Soils
At our study site, we observed that sites closer to the stream tended to have higher P concentrations than at sites at greater distance to it (Figures 1 and 2 Table S4), and the green shading represent the linear mixed effects model 95% confidence interval. The gray band indicates the range of equilibrium values for isotopic exchange between phosphate and water, which was calculated according to Chang and Blake (2015) using water isotope and temperature data for the site and accounting for a ±2‰ uncertainty in water isotope values due to evaporative enrichment of the soil water (see Methods section for details) et al., 2009). In addition, all sites were characterised by high concentrations of P o , supporting recent suggestions that P o constitutes a major proportion of P transported through the VBS (Badiou et al., 2018;Dodd et al., 2018;Young et al., 2013). Our study site was managed as grassland in the last few decades, and no excessive flooding has been observed during the past 25 yr. Hence, we discount the possibility that P patterns were influenced by recent flooding and sedimentation events. Instead, the elevated P concentrations in soils near the stream are likely a legacy of deposition of material eroded from the surrounding hillslopes, which is supported by the colluvial properties (AB horizons) observed in the soil profiles. We therefore interpret the trend of higher P concentrations at locations closer to the stream as an indication of long-term retention of soil and nutrients within the VBS (S. Gu et al., 2017;Neidhardt et al., 2019;W. Roberts et al., 2012). From the limited number of studies that considered different soil depths, it has been suggested that P accumulation in VBS soils only occurs at the soil surface (Aye et al., 2006;S. Gu et al., 2017). In contrast to these previous findings, our data showed high P concentrations also in the subsoil of the VBS. Also, including sampling depth improved our ability to predict P concentrations (Supplemental Tables S2 and S3). On the one hand, high P concentrations at the soil surface may reflect deposition of eroded material over several cen-turies, which may now serve as additional sources of P release to freshwater bodies (Weihrauch & Weber, 2020;Weihrauch et al., 2021). On the other hand, deeper soil horizons have been identified as potential sinks for P (Aye et al., 2006), indicating that P retained from overland flow or erosion may have been partially transferred to deeper soil horizons by leaching or particulate transport. Our data point to a need to consider deeper soil horizons in a characterization of P retention and cycling in VBS systems. It should be noted that linear models including the interaction term of distance and depth provided a better fit for our data than a simple additive model (Supplemental Tables S2 and S3). This suggests that the tendency for P accumulation along the transect from grassland to VBS may have been influenced by additional nonadditive effects of environmental parameters, such as P discharge, annual or seasonal water levels and redox potentials, which were not recorded in this pilot study.
Any P discharge would either involve P bound to soil particles (Gottselig et al., 2020;Young et al., 2021) or labile P fractions such as H 2 O-and 0.5 M NaHCO 3 -extractable P. Indeed, the labile P fractions comprised only a minor fraction of extracted P along the transect, suggesting that at least some of this labile P might have been transferred into the stream, despite the presence of the VBS. The 0.1 M NaOH P fraction appears to be affected most by the changing soil moisture conditions along the transect towards the river as it tended to have up to threefold higher P concentration at the site close to the river compared with site furthest away. The 0.1 M NaOH P fraction extracts less available P bound to Fe or Al oxides (Cross & Schlesinger, 1995;Tiessen & Moir, 1993). Close to the stream, the seasonal fluctuation of the groundwater table and associated capillary rise results in frequent changes of oxidation and reduction processes (S. Gu et al., 2017). Changing redox conditions induce a shift in Fe chemistry, with insoluble Fe 3+ dominating under aerobic conditions and easily soluble Fe 2+ under anaerobic conditions (S. Gu et al., 2017;Reddy et al., 1999). This dissolution and precipitation of Fe-oxides also releases P into or precipitates P from the soil solution (Parsons et al., 2017). Thus, the 0.1 M NaOH-extractable P pool should be the main fraction affecting P source or sink functions in VBS soils.
The P concentration in the 1 M HCl fraction changed only marginally along the transect (Figures 1 and 2); distance alone described <5% of the variance in P i and P o concentrations with depth (Supplemental Tables S2 and S3). There is evidence that Fe from poorly crystalline oxides can also be extracted by 1 M HCl (C. Gu et al., 2020;Slomp et al., 1996), which may also transfer Fe-associated P into the HCl extract. More often it is assumed that 1 M HCl-extractable P includes P in Ca-minerals, either in the form of primary minerals inherited from the parent material or secondary minerals newly precipitated from free Ca-and phosphate-ions (Angert et al., 2011;Helfenstein et al., 2018).
In summary, the trends that we observed for the distribution of the soil P pools along the VBS to grassland transect agree with previous observations on VBS soils across climates and ecosystems. However, although previous studies have also suggested that water dynamics and P accumulation should affect biological P cycling in VBS soils (Stutter et al., 2009;Young et al., 2021), direct evidence has been lacking so far.

Biological cycling of P in VBS soils
Oxygen isotopes in phosphate have been extensively used to study biological P cycling in soils (e.g., Helfenstein et al., 2018;Joshi et al., 2016;Tamburini et al., 2012) and for source attribution of P in water bodies (Elsbury et al., 2009;Ishida et al., 2019;Paytan et al., 2017;Pistocchi et al., 2017;Tonderski et al., 2017). The oxygen isotope values in the soil profiles along our transect generally fall into the range of isotope values observed in other studies in temperate climates (Amelung et al., 2015;Bauke et al., 2018;K. Roberts et al., 2015). However, although 1 M HCl-extractable P concentrations changed only marginally along the transect (Figures 1  and 2), we observed a trend of decreasing δ 18 Ο P values with increasing distance from the stream (Figure 3). We attribute this finding to higher intensity of biological soil P cycling close to the stream than in the adjacent grassland soils. The δ 18 Ο P values may thus provide a more sensitive tracer of processes occurring in VBS soils than singular or repeated analyses of P fraction concentrations. The observation of soil δ 18 Ο P values closer to equilibrium near the stream compared with values from soils that were located further away, is in line with earlier reports of elevated microbial activity and respiration in VBS compared with adjacent agricultural soils (Jacinthe et al., 2015;Stutter et al., 2009;Tufekcioglu et al., 2001). Thus, we suggest that our findings provide indications of increased microbial P cycling in VBS soils and soils that undergo fluctuating redox conditions (Bernhardt et al., 2017;Saia et al., 2021).
In topsoils, the δ 18 Ο P value was substantially higher than in deeper parts of the soil profiles. Apart from more intense biological P cycling in A than in C horizons, elevated δ 18 Ο P values could also derive from plant litter that has accumulated in the topsoil and which is typically enriched in 18 O (Bauke et al., 2021;Pfahler et al., 2013Pfahler et al., , 2017. Especially for the site 1 m away from the stream, frequent plant litter inputs from the dense VBS vegetation and subsequent hydrolysis of associated P o compounds by the activity of extracellular acid phosphatases provide a source of 18 O-enriched phosphate, which can then be transferred to the P i pools in the soil (Helfenstein et al., 2018;von Sperber et al, 2014).
The additive effect of distance to the stream and soil depth represented the best model explanation for the varia-tion in δ 18 O P values (Supplemental Table S4). Our findings thus agree with previous suggestions indicating that nutrient cycling is not limited to the topsoil and may extend into the subsoil (Amelung et al., 2015;Bauke et al., 2018). With greater distance to the stream, the δ 18 O P values decreased first in the subsoil and eventually also in the topsoil at the most remote site (Figure 3), pointing to a trend of decreasing biological P cycling. Considering that the subsoil of the VBS and grassland close to the stream are most affected by fluctuating water tables, the reductive dissolution of Fe-bound phosphate and the biological cycling of P in these subsoils may present an additional source of P mobilization and release to the aquifer that has so far been neglected (Weihrauch & Weber, 2020;Weihrauch et al., 2021). Future studies should therefore aim to quantify the spatial extent of the riparian zone not only across the landscape, but also with respect to the depth of P storage and cycling as well as soil redox conditions.
As a first estimate of the influence of redox conditions on biological P cycling, we classified our data into those belonging to the topsoil Ah (i.e., permanently oxidized) and the subsoil Cr horizon (i.e., permanently reduced). We found that the additive effect of transect distance and redox condition explained 82% of variation in δ 18 Ο P values (Supplemental Figure S3 and Table S5). Given the limited number of samples included here, additional research may be needed to test to what extent differences in microbial metabolism in oxic and anoxic conditions affect oxygen isotope exchange processes. Validating our findings for other riparian zones therefore now warrants further attention.
This study was, to the best of our knowledge, the first that analyzed δ 18 Ο P values in VBS soils. The δ 18 Ο P values of the 1 M HCl pool can be assumed to integrate information on P cycling over time scales of months to years (Angert et al., 2012;Helfenstein et al., 2018;Joshi et al., 2016). Over these time scales, a small but continuous fraction of P cycled in organismic cells may be precipitated into secondary minerals, which will slowly imprint the signal of microbial cycling to HCl-extractable P. Based on the trends in P cycling observed here even in this stable P pool, follow-up studies should be carried out with a higher spatial and temporal resolution and consider more dynamic P pools to obtain a better understanding of the extent and mechanisms of P cycling in VBS soils.
As a consequence of enhanced P cycling in VBS soils, it is unlikely that the P released into the stream carries the isotopic value of the source material. Rather it is more likely that the P released into the stream carries the isotopic equilibrium value that occurred in VBS soils. So far, studies that use oxygen isotope data to identify P sources to surface waters (Elsbury et al., 2009;Ishida et al., 2019;Paytan et al., 2017;Pistocchi et al., 2017;Tonderski et al., 2017) have done so by only comparing the original source values of soil and fertilizer with those observed for phosphate in water but hardly considered the biological cycling of P in VBS soils or in proximity of the river banks. We suggest that a quantitative assessment of P source contributions can only be achieved if biological transformations along the transport path are considered, including P cycling in VBS soils.

CONCLUSION
The accumulation of P in the temperate VBS soils studied here was linked to enhanced biological P cycling. Specifically, we observed a trend of enhanced biological P cycling near the stream, which modified the oxygen isotopic composition of phosphates compared with the adjacent grassland soils. Hence, the results of this study indicate that P transformation along the transport pathway from soil to stream may affect P source attribution and budgeting estimates. This effect is currently not considered in studies using oxygen isotopes as a tracer of P inputs into surface waters. Elucidating the factors involved in this biological P cycling in VBS soils, such as pore size, water chemistry, soil redox conditions, P saturation status, and microbial P limitation, now warrants further attention.

ACKNOWLDGMENTS
We would like to thank Dr. Jens Kruse for helping with the analysis of phosphorus concentrations in extracts, Dr. Patrick Sullivan for help with statistical analyses, and two anonymous reviewers for their constructive feedback.

D A T A AVA I L A B I L I T Y S T A T E M E N T
All data and analysis scripts associated with this publication are available on GitHub (https://github.com/sheilasaia/paper-d18o-phosphate-riparian) and Zenodo (DOI: 10.5281/zenodo.5920850, https://zenodo.org/record/5920851#.Yfb04mB-OlYs).

C O N F L I C T O F I N T E R E S T
The authors declare no conflict of interest.