Introduction

In recent years, the potential environmental impact of high concentrations of inorganic phosphate (PO4) in groundwater ecosystems has received increasing attention1,2,3. For instance, groundwater with high PO4 concentrations can be transported to surface waters and induce eutrophication4,5. Moreover, high concentrations of dissolved PO4 in groundwater can promote the mobility of toxic arsenic via competitive adsorption6,7, which poses a threat to groundwater quality and human health, especially in SE Asia8. Therefore, it is important to reveal the biogeochemical processes that are involved in the release of dissolved PO4 in groundwater.

Previous studies focused on the anthropogenic inputs of PO4 into groundwater at contaminated sites (i.e., farmlands9,10,11 and wastewater treatment plants12,13,14), but recent studies highlighted the importance of geogenic, aquifer-internal PO4 release in groundwater of floodplain aquifers worldwide15,16,17. Generally, organic matter (OM) degradation, reductive dissolution of Fe(III) oxides and apatite dissolution are regarded as the three main processes that release PO4 into the groundwater18,19,20. In aquatic ecosystems, other biological processes may also contribute to PO4 release. Under natural conditions, enzymatically mediated reactions release PO4 from organic phosphorus compounds (Porg) to meet the demand of microbial growth21. The regenerated PO4 eventually undergoes microbial uptake and intracellular cycling, and is finally released back to the extracellular environment21. Previous studies demonstrated that microbial cycling is an important process that releases PO4 in soil, lakes, and marine environments22,23,24,25,26. In floodplain aquifers, it is known that microbial activity is linked to OM degradation and Fe(III) reduction27. However, whether released PO4 has been utilized by microorganisms and how it contributes to the increase of dissolved PO4 in floodplain aquifers remains unclear.

In soils and marine environments, it is known that microbial activity and PO4 utilization are affected by changing redox conditions28,29,30. On the one hand, redox conditions directly constrain the rate of microbial respiration and activity, therefore the ability of microbes to utilize PO429,31. On the other hand, microbial PO4 utilization is indirectly influenced by redox conditions30,32. Under anoxic conditions, the reductive dissolution of Fe(III) oxides can release surface-adsorbed, occluded and/or co-precipitated PO433,34, thus potentially enhancing the PO4 bioavailability and microbial PO4 utilization. Previous studies in floodplains have demonstrated tight constraints of redox conditions on microbial activity35,36 and PO4 release in aquifers37. Furthermore, Fe(III) oxides are abundant in aquifer sediments serving as an important sink for PO47, which might eventually be released under anoxic conditions. Therefore, the contribution of microbially cycled PO4 could either be decreased or increased under anoxic as compared to oxic conditions in floodplain aquifers.

In the past decade, the stable isotope ratio of PO4-O (δ18OPO4) has emerged as a useful tool for identifying PO4 sources38,39,40,41 and the involvement of microbial activity in PO4 cycling in surface aquatic ecosystems42,43,44. Under typical groundwater temperature, pressure, and pH conditions, the P-O bonds in PO4 are resistant to inorganic hydrolysis45,46. However, enzymatic cleavage of PO4 is associated with an exchange of O atoms from ambient water accompanied by isotopic fractionation47,48,49. In particular, the intracellular enzyme pyrophosphatase imposes a temperature-dependent equilibrium isotopic fractionation of O between PO4 and ambient water50,51. As a result, δ18OPO4 values may shift from the original source signature (i.e., minerals in the bedrock from which the PO4 was originally released) towards isotopic equilibrium values (δ18OPO4 equ)42,52,53. Theoretical δ18OPO4 equ values are calculated for a given sample based on the local temperature and δ18O of ambient (ground)water (δ18OH2O)50. So far, only a few studies have investigated the δ18OPO4 characteristics of groundwater54 and their results showed that δ18OPO4 values may differ from the equilibrium value, possibly reflecting anthropogenic PO4 inputs15,55. However, in floodplain aquifers, high concentrations of PO4 and other typical geogenic contaminants are most likely accumulated from long term water-rock interactions20,56,57. Therefore, δ18OPO4 characteristics in groundwater of floodplain aquifers may behave differently, reflecting aquifer-internal PO4 release processes.

The Hetao Basin is located in the Inner Mongolia Autonomous Region, PR China, where more than 70% of the groundwater has phosphorus (P) concentrations exceeding 50 μg L−1, and dissolved P in groundwater predominantly presents in the form of PO437. In the aquifer sediments, Fe(III) oxides represent an important source of PO4 with average contents of adsorbed PO4-P of 9.1 mg kg−1. The reductive dissolution of Fe(III) oxides is considered as a dominant process releasing PO4 under anoxic conditions37. Along the groundwater flow path, redox conditions gradually change from oxic to strongly reducing which is associated with a sequence of microbially mediated redox reactions, such as microbial iron oxidation, microbial Fe(III) reduction and microbial sulfate reduction58,59,60 (Fig. 1). There are multiple clay layers in the sediments, separating a shallow (5–40 m bls) and a deep aquifer (40–100 m bls)61. Our previous investigation showed that dissolved PO4 concentrations were significantly higher in the deep aquifer than in the shallow aquifer, suggesting a dominant role of in situ geochemical processes in the release of PO437. Moreover, we found that concentrations of NH4+, NO3, and SO42−, being potential proxies for anthropogenic inputs10,12, were all significantly lower in the shallow aquifer than in the deep aquifer37. Hence, in this study, we consider a separation of the aquifers from the surface and apply the δ18OPO4 analysis to (1) uncover microbial cycling of PO4 in groundwater systems and (2) reveal the influence of redox conditions on microbial PO4 cycling in groundwater of floodplain aquifers.

Fig. 1: Groundwater sampling locations within the north-western Hetao Basin, Inner Mongolia Autonomous Region, PR China.
figure 1

The dark blue arrow shows the approximate groundwater flow direction from the Langshan Mountains in the north-west towards the flat plain in the centre of the Hetao Basin. Along the flow path, 13, 16, and 16 groundwater samples were collected in the alluvial fan (the blue circles), the transition area (the green circles) and the flat plain (the red circles), respectively. In particular, the bold circles show sampling locations of 12 groundwater samples that have been previously analyzed for 16S rRNA gene amplicon sequencing60,113,114, including four samples in the alluvial fan (LY07, LY08, 17-7-11, and K2-6), four samples in the transition area (I-9, H01-30, I-3, and I-1), and four samples in the flat plain (K1-2, K1-6, I-16, and 1803). Maps created in Esri ArcMap 10.5 (ESRI World Image Basemap, source: Esri, Maxar, Earthstar Geographics).

We found that the contribution of microbial PO4 cycling to the dissolved PO4 significantly increased along a groundwater flow path from oxic to anoxic conditions. In the alluvial fan under oxic conditions, although 25–47% of the dissolved PO4 inherited the initial source signal of igneous apatite, groundwater δ18OPO4 reflected a pronounced impact of intracellular enzymatic cycling by microorganisms. In the flat plain under anoxic conditions, dissolved PO4 carried a nearly exclusive equilibrium isotope signal. This was probably due to the release of PO4 with an equilibrium δ18OPO4 value from Fe(III) oxides, as a result of Fe(III) reduction in the presence of Fe(III)- and sulfate-reducing bacteria. Another potentially important aspect was the increasing residence times of groundwater from the alluvial fan ( < 125 years) to the flat plain (ca. 600 years), which allowed the dissolved PO4 to be more microbially cycled and consequently the groundwater δ18OPO4 values gradually approached isotopic equilibrium. We therefore conclude that PO4 in groundwater is tightly microbially cycled under a wide range of redox conditions, and highlight microbial PO4 cycling involved in the reductive dissolution of Fe(III) oxides as an overlooked biological process contributing to the release of PO4 into groundwater.

Results and discussion

The presence of microbial PO4 cycling in groundwater

Along the groundwater flow path, the δ18OH2O of groundwater increased from the alluvial fan (median of −10.6‰) over the transition area (median of −9.9‰) towards the flat plain (median of −8.7‰) (Figs. 2a, 3a). In the basin, groundwater levels gradually decreased from the alluvial fan (about 10 m bls) to the flat plain (about 2 m bls)56. The gradually increased δ18OH2O values most likely reflected increasing evaporation along the flow path. This was further supported by the increasing TDS values and Cl concentrations (see Supplementary Table S1), agreeing with previous studies58,62. Similar to δ18OH2O, δ18OPO4 values in groundwater significantly increased (p < 0.01) from the alluvial fan (median of +11.7‰) over the transition area (median of +13.4‰) towards the flat plain (median of +15.2‰) (Figs. 2b, 3b). Since abiotic processes have negligible effects on the δ18OPO4 in aquatic ecosystems63,64, the increase in δ18OPO4 values along the flow path and their significant positive correlation with δ18OH2O (R2 = 0.44; p < 0.01) (Fig. 4) suggest a pronounced influence of intracellular enzymatic cycling of P, resulting in an exchange of O atoms between PO4 and ambient groundwater.

Fig. 2: Spatial variation in groundwater δ18OPO4 values and redox-sensitive components along the flow path.
figure 2

a Observed δ18O isotope values of groundwater (δ18OH2O). b Observed δ18O isotope values of dissolved inorganic phosphate in groundwater (δ18OPO4). c Calculated theoretical equilibrium δ18OPO4 values of dissolved inorganic phosphate after intracellular microbial cycling (δ18OPO4 equ). d Difference between observed δ18OPO4 and δ18OPO4 equ values (Δ18OPO4) (Negative Δ18OPO4 values indicate that observed δ18OPO4 were isotopically lighter (i.e. more 18O-depleted) than the calculated equilibrium δ18OPO4 equ values). e Redox potential (Eh) of groundwater. fh Concentrations of dissolved inorganic phosphate (PO4), ferrous iron (Fe(II)), and sulfide (S2-), respectively. The blue, green and red rectangles showed the ranges of the alluvial fan, the transition area and the flat plain, respectively. Results of statistical analyses are shown in Fig. 3.

Fig. 3: Statistical analyses on groundwater δ18OPO4 values and redox-sensitive components in the different zones.
figure 3

a Observed δ18O isotope values of groundwater (δ18OH2O). b Observed δ18O isotope values of dissolved inorganic phosphate in groundwater (δ18OPO4). c Calculated theoretical equilibrium δ18OPO4 values of dissolved inorganic phosphate after intensive intracellular microbial cycling (δ18OPO4 equ). d Difference between δ18OPO4 and δ18OPO4 equ values (Δ18OPO4) (Negative Δ18OPO4 values indicate that observed δ18OPO4 were isotopically lighter (i.e. more 18O-depleted) than the calculated equilibrium δ18OPO4 equ values). e Redox potential (Eh) of groundwater. fh Concentrations of dissolved inorganic phosphate (PO4), ferrous iron (Fe(II)), and sulfide (S2-), respectively. The Zone I, II and III stand for the alluvial fan, the transition area and the flat plain, and are shown in blue, green, and red colours, respectively. Uppercase letters (A–C) indicate significant differences (p < 0.05) in measured variables among the three zones. N = 13, 16, and 16 in the alluvial fan, the transition area, and the flat plain, respectively.

Fig. 4: Relationship between observed δ18OPO4 and δ18OH2O values in groundwater.
figure 4

δ18OPO4 values represent observed δ18O isotope values of dissolved inorganic phosphate in groundwater; δ18OH2O values represent observed δ18O isotope values of groundwater. The regression coefficients (R2) are provided together with the level of significance (p).

In natural ecosystems, the intracellular enzymatic cycling of PO4 in microorganisms is associated with a temperature-dependent equilibrium isotope fractionation65. In the study area, the groundwater temperatures showed only a small variation (10.9 ± 0.11 °C, see Supplementary Table S1). Based on Eq. (1) (see Data processing in Methods), we calculated the theoretical δ18OPO4 equ values that would be expected in the case of a complete overprinting by intracellular microbial PO4 cycling50. Our results showed that δ18OPO4 equ values significantly increased along the flow path (medians of +13.6‰, +14.4‰ and +15.6‰ in the alluvial fan, the transition area and the flat plain, respectively; Figs. 2c, 3c). More importantly, the difference between observed δ18OPO4 values and expected δ18OPO4 equ values (Δ18OPO4, see Eq. (2) in Methods) was relatively small, with a median value of less than −1‰ in the studied groundwater (Figs. 2d, 3d). This small difference suggests that O isotope ratios of dissolved PO4 molecules in groundwater were significantly affected by microorganism intracellular cycling. Similarly, microbial PO4 cycling was observed in other aquatic ecosystems, such as marine sediments and lakes, as indicated by δ18OPO4 values close to δ18OPO4 equ values52,66. Hence, microbial cycling of PO4 also plays an important role in groundwater systems.

Increasing contribution of microbial cycling to dissolved PO4 from oxic to anoxic conditions

Along the groundwater flow path, the Eh values sharply decreased from the alluvial fan (median of 137 mV) to the transition area (69.4 mV) and the flat plain (60.9 mV; Figs. 2e, 3e), indicating that the redox conditions of groundwater changed from oxic to anoxic. The median concentrations of dissolved PO4 significantly (p < 0.05) increased from the alluvial fan (51.0 μg L−1) to the transition area (109 μg L−1) and the flat plain (92.5 μg L−1; Figs. 2f, 3f), showing an increase of dissolved PO4 concentrations under anoxic conditions. Furthermore, Δ18OPO4 values in groundwater significantly increased (p < 0.01) from the alluvial fan (median of −2.0‰), over the transition area (median of −1.1‰), to the flat plain (median of −0.2‰), approaching zero (Figs. 2d, 3d). This indicated that the difference between observed δ18OPO4 values and expected δ18OPO4 equ values was getting smaller along the flow path. The significant negative correlation between Δ18OPO4 and Eh values (r = −0.50; p < 0.01; Fig. 5a) indicated an increasing contribution of microbial PO4 cycling to the dissolved PO4 along with the change from oxic to anoxic conditions. Hence, biogeochemical processes that are potentially associated with microbial PO4 cycling and result in the changes of Δ18OPO4 values from oxic to anoxic conditions are discussed as following.

Fig. 5: Relationship between Δ18OPO4 values and redox-sensitive components along the groundwater flow path.
figure 5

Δ18OPO4 values represent the isotope difference between observed δ18O isotope values of dissolved inorganic phosphate in groundwater (δ18OPO4) and calculated theoretical equilibrium δ18OPO4 values of dissolved inorganic phosphate after intracellular microbial cycling (δ18OPO4 equ). ad Δ18OPO4 values against redox potential values (Eh), concentrations of dissolved ferrous iron (Fe(II)), dissolved inorganic phosphate (PO4), and dissolved sulfide (S2-) in groundwater, respectively. The blue, green and red circles represent groundwater samples in the alluvial fan, the transition area and the flat plain, respectively. The correlation coefficients (r) are provided together with the level of significance (p).

Microbial PO4 cycling under oxic conditions

Although microbial PO4 cycling contributed to the dissolved PO4, the largest Δ18OPO4 (median value of −2.0‰; Fig. 3d) and the most 18O-depleted isotope ratios of PO4 in the alluvial fan suggest that δ18O of dissolved PO4 in groundwater under oxic conditions was closer to its source signal. Based on our previous investigation37, either sedimentary organic matter (SOM) and/or apatite were considered as potential sources of PO4 under oxic conditions in our groundwater system.

In the study area, SOM in aquifers is degradable67 and mostly originates from buried plant residues in the unconsolidated sediments68. Typical δ18OPO4 values of plants range from +17‰ to +27‰69, which is reasonable to expect for the SOM as well. Previous experiments demonstrated that the extracellular Porg hydrolysis commonly occurs in two steps, which are catalyzed by phosphodiesterases (diesterases) and phosphomonoesterases (monoesterases), respectively47. For each step, PO4 inherits three O atoms from the Porg molecule and incorporates one O atom from ambient water70,71. Hence, PO4 released from extracellular Porg hydrolysis contains signatures of both enzymatic steps and the combined isotopic effect of diesterase + monoesterase is expressed as Eq. (3)47 (see Data processing in Methods). Considering the weakly alkaline pH of groundwater in the alluvial fan (pHs: 6.7 – 9.1; see Supplementary Table S1), alkaline phosphatase (APase) can be assumed to be the key extracellular monoesterase enzyme that was involved in OM degradation and the concomitant release of PO4 to groundwater49,72. In natural waters, PO4 is expected to have a stronger RNA degradation signature than DNA degradation signature47. For RNA, the fractionation factor between O incorporated into PO4 and ambient water O is −5( ± 6)‰ for the diesterase + APase pathway70,71. Therefore, δ18OPO4 values of PO4 from SOM biodegradation catalyzed by extracellular diesterase and APase enzymes would range from +0.0‰ to +5.5‰, shifting the observed δ18OPO4 values away from the δ18OPO4 equ values towards more 18O-depleted signatures22,69,73. However, NH4+ concentrations in groundwater of the alluvial fan were the lowest (p < 0.05) among the three zones with median value of 0.40 mg L−1 (see Supplementary Table S1). In groundwater ecosystems, since NH4+ mainly originates from the mineralization of N-containing OM5, NH4+ represented a useful indicator of OM degradation74,75. According to the lowest concentration of NH4+, OM degradation played a less important role in PO4 release in the alluvial fan as compared to the other two zones. This was in agreement with our previous investigations, demonstrating that when compared with the transition area and the flat plain, the bioreactivity of OM and the Porg contents in aquifer sediments of the alluvial fan were relatively lower35,37. Hence, PO4 released from OM degradation might not result in the largest Δ18OPO4 observed in the alluvial fan. Therefore, we infer a limited contribution of OM degradation to PO4 concentrations in groundwater under oxic conditions.

In the aquifers, apatite dissolution is considered another source of dissolved PO437,76. In the Hetao Basin, we found that Ca/Fe(II)-minerals represented the dominant sedimentary PO4 pool and also that groundwater in the alluvial fan was undersaturated with regard to apatite37. We further measured δ18OPO4 values of HCl-PO4 (representing PO4 in apatite minerals) in the aquifer sediments, which ranged from +4.0‰ to +8.5‰ (see Supplementary Table S2). The measured δ18O values of HCl-PO4 were in line with those previously reported for soils and sediments derived from igneous rocks (from +5.3‰ to +9.0‰)77,78. In the adjacent Langshan Mountains, the contents of PO4-P in metamorphic complex were reported to exceed 1000 mg kg−179. Notably, there were two groundwater samples with remarkably low δ18OPO4 values (+5.9‰ and +8.7‰) in the alluvial fan (Fig. 2b) that matched δ18OPO4 values of HCl-PO4 in the sediments (from +4.0‰ to +8.5‰). This indicated that the release of PO4 from igneous apatite dissolution characterized by low δ18OPO4 signals was the main process shifting δ18OPO4 values of groundwater away from isotopic equilibrium.

In the alluvial fan, the mean of observed δ18OPO4 values in groundwater was +11.2‰ (see Supplementary Table S1). Based on a simple two end-member mixing model with the measured δ18OPO4 of igneous apatite (with +4.0‰ as lower and +8.5‰ as upper boundary) and the calculated δ18OPO4 equ of groundwater (mean value of +13.6‰ in the alluvial fan, see Supplementary Table S1), we estimated that 25–47% of the dissolved PO4 was derived from igneous apatite dissolution. We therefore conclude that apatite dissolution contributed to the isotopic signal of PO4 in groundwater under oxic conditions, whereas > 53% of the dissolved PO4 was intensively turned over by microorganisms. Previous studies showed that in soil and marine environments, microbial activity and turnover releases PO4 into solution under oxic conditions23,25,80,81. Thus, our results highlight that dissolved PO4 in oxic aquifers may be subject to active microbial cycling of PO4.

Microbial PO4 cycling under anoxic conditions

As previously described, dissolved PO4 was increased under anoxic conditions in the transition area and the flat plain. Along the groundwater flow path, Fe(II) concentrations significantly (p < 0.01) increased from the alluvial fan (0.33 mg L−1 in median) over the transition area (0.84 mg L−1 in median) to the flat plain (1.94 mg L−1 in median) (Figs. 2g, 3g). Since Fe(III) oxides represent a potential source for PO4 in the aquifer sediments from the basin (average: 9.10 mg kg−1)37, the positive correlation between groundwater Fe(II) and PO4 (r = 0.53, p < 0.01) suggests that the reductive dissolution of Fe(III) oxides was an important process of PO4 mobilization and resulted in the high concentrations of dissolved PO4 under anoxic conditions. Simultaneously, the Δ18OPO4 values gradually increased along the groundwater flow path, approaching zero (Fig. 3d). Especially in the flat plain, the median Δ18OPO4 value was not significantly different from zero (p > 0.05), reflecting that O isotope equilibration was reached between the dissolved PO4 and ambient groundwater. Hence, under anoxic conditions, it appears that the isotopic signal of the dissolved PO4 mostly originates from microbial cycling. There are two possible explanations for the observed changes of Δ18OPO4 value along the flow path. One is that more PO4 with an equilibrium δ18OPO4 value was released into groundwater via the reductive dissolution of Fe(III) oxides under anoxic conditions, and the other is that dissolved PO4 was more intensively cycled by microorganisms in groundwater of the flat plain.

Regarding the first explanation, it would require that PO4 associated with Fe(III) oxides attained an isotopic equilibrium value before and/or during the reductive dissolution of Fe(III) oxides. We cannot exactly determine the time when the microorganism left their isotope imprint on Fe(III)-associated PO4. Some studies showed that in a thermally stratified Archaean ocean, Fe(III)-associated PO4 recorded equilibrium δ18OPO4 values with cooler ocean surface waters52. Moreover, studies on Fe(III) oxide deposits from both Red and Green seamounts located near the East Pacific Rise (EPR) concluded that PO4 bound to Fe(III) oxides was intensively cycled by Fe-oxidizing bacteria as evidenced by mineral morphology and δ18OPO4 values of PO4 in Fe(III) oxides54. Specifically, Fe(III) oxides in both seamounts were mainly composed of twisted filaments, stalks, and hollow sheaths that resembled the characteristic remains of the Fe-oxidizing bacteria Gallionella ferruginea and Leptothrix82. In the study area, Fe-oxidizing bacteria were identified and dominated in groundwater of the alluvial fan with an averaged relative abundance of 7.6%, including Gallionella83, Acidovorax84, Thermomonas85, an unclassified member of the Gallionellacae family and an unclassified member of the Rhodobacteraceae86 family (Fig. 6). Moreover, the enriched functional genes of foxEY, Cyc1 and sulfocyanin were observed in the alluvial fan, further indicating that microbial Fe oxidation was the major Fe-related biogeochemical process87. Hence, microbial imprint might have happened before the reductive dissolution of Fe(III) oxides, and therefore PO4 associated with Fe(III) oxides attained an equilibrium value under oxic conditions. In this case, if PO4 was later released via the reductive dissolution of Fe(III) oxides without further microbial imprint, we expect that δ18OPO4 values in groundwater under anoxic conditions were similar to δ18OPO4 equ values under oxic conditions. However, as we observed, δ18OPO4 values in the flat plain under anoxic conditions ( + 14.0‰ - +16.7‰) were significantly higher than δ18OPO4 equ values in the alluvial fan under oxic conditions ( + 13.1‰ - +14.0‰) (see Supplementary Table S1), which suggest that microbial imprint on PO4 associated with Fe(III) oxides was further generated under anoxic conditions. This was consistent with the findings that both dissolved PO4 in porewater and Fe(III) oxide-associated PO4 in anoxic bottom sediments in the Chesapeake Bay had δ18OPO4 values near equilibrium values88, and that the release of PO4 from Fe(III)/Al oxides in lake sediments under anaerobic conditions led to high PO4 concentrations in Lake Shijiuhu in the Yangtze River Basin, with biotic activity regulating the δ18OPO4 values toward equilibrium41.

Fig. 6: Microbial community composition along the groundwater flow path.
figure 6

From the alluvial fan to the flat plain, relative abundance of iron-oxidizing bacteria (bars in blue series) gradually decreased, while relative abundances of Fe(III)-reducing bacteria (bars in red series) and sulfate-reducing bacteria (bars in yellow series) gradually increased, showing an shift of microbial community from oxic to anoxic conditions.

When redox conditions gradually changed from oxic to anoxic along the groundwater flow path, Fe(III)-reducing bacteria became the dominant microorganism that utilized PO4 associated with Fe(III) oxides. This was evidenced by the higher abundance of Fe(III)-reducing bacteria in the flat plain (Fig. 6), including Aeromonas89, Geobacter90, Paludibaculum91, Rhodoferax92 and Sphingomonas93. Along the flow path, the relative abundance of those dominant Fe(III)-reducing bacteria remarkably increased from the alluvial fan (0.7% in average), over the transition area (2.3% in average), to the flat plain (2.9% in average) (Fig. 6). Furthermore, Xiu et al.87 found that the enriched functional gene of omcS was present in the same flat plain with the abundance of 0.6 fragments per kilobase per million mapped fragments (FPKM), indicating microbial Fe(III) reduction became the major Fe-cycling process. Previous experiments showed that a significant fraction of Fe(III) oxides-bound PO4 was bioavailable and was taken up by bacteria attaining an isotopic equilibrium composition53. Moreover, it is well-known that Fe(III)-reducing bacteria can access solid Fe(III) as electron acceptor, inducing the reductive dissolution of Fe(III) oxides via multiple electron transfer strategies94. As we observed, along the groundwater flow path, Fe(II) concentrations and Δ18OPO4 value consistently increased with a significantly positive correlation (r = 0.47, p < 0.01; Fig. 5b). Simultaneously, the concentrations of PO4 also increased with Δ18OPO4 value with a positive correlation (r = 0.34, p < 0.05; Fig. 5c). The significant positive correlations between Δ18OPO4 value and dissolved Fe(II), and PO4 suggest that more PO4 with an equilibrium δ18OPO4 value was released into groundwater via the reductive dissolution of Fe(III) oxides, as a result of increasing microbial Fe(III) reduction under anoxic conditions in the flat plain (Fig. 7). Similarly, Zhao et al. observed that a shift of δ18OPO4 towards equilibrium at a depth around 140 m below the seafloor in the Peru Margin upwelling zone corresponded well with a shift in microbial communities44.

Fig. 7: Conceptual figure of the increasing contribution of microbial cycling to dissolved inorganic phosphate (PO4) from oxic to anoxic conditions in floodplain aquifers of the Hetao Basin.
figure 7

The green points show the observed δ18O isotope value of PO418OPO4) in groundwater. The red rectangle stands for the range of equilibrium δ18OPO4 values after intracellular microbial cycling. The blue rectangle represents the range of δ18OPO4 values of igneous apatite. Along the flow path, the redox conditions gradually changed from oxic to anoxic, and the groundwater residence time increased from about 125 yr to 675 yr. Under these conditions, the green dashed line shows an increasing trend of δ18OPO4 values towards equilibrium δ18OPO4 values along the groundwater flow path, indicating an increasing contribution of microbial cycling to dissolved PO4 from 53% − 75% to nearly 100%. This is probably due to the release of PO4 with an equilibrium δ18OPO4 value from Fe(III) oxides, as a result of Fe(III) reduction in the presence of Fe(III)- and sulfate-reducing bacteria (red dashed lines). Moreover, the increasing residence times of groundwater allow the dissolved PO4 to be more microbially cycled and consequently the groundwater δ18OPO4 values gradually approached isotopic equilibrium.

Under anoxic conditions, microbial sulfate reduction would further promote the reductive dissolution of Fe(III) oxides with attendant release of PO4 with an equilibrium δ18OPO4 value into groundwater. Similar to Fe(III)-reducing bacteria, sulfate-reducing bacteria also utilized PO4 associated with Fe(III) oxides under anoxic conditions. This is evidenced by the presence of available sulfate (median of 225 mg L−1), and sulfate-reducing bacteria (Fig. 6) (including Desulfovibrio95, Desulfuromonas96, Desulfatirhabdium97, an unclassified member of the Desulfosarcinaceae95 family, and an unclassified member of the Desulfobacteraceae95 family) in groundwater of the transition area and the flat plain. Specifically, the averaged relative abundance of those dominant sulfate-reducing bacteria increased from the alluvial fan (0.2%), over the transition area (0.6%), to the flat plain (1.3%) (Fig. 6). In addition, in our previous study87, the enriched functional genes related to dissimilatory sulfate reduction (average FPKM from 0.4 to 1.1), such as sat, APR, aprAB, and dsrAB, were observed in groundwater of the same transition area and flat plain. As shown in Figs. 2h and 3h, the concentrations of S2- in the transition area (median of 50.0 μg L−1) were significantly higher than those in the alluvial fan (median of 25.0 μg L−1), which indicated that microbial sulfate reduction occurred in groundwater under anoxic conditions. Both S2- and Fe(II) increased from the alluvial fan to the transition area with a positive correlation (r = 0.57, p < 0.01), suggesting that the reductive dissolution of Fe(III) oxides was promoted by microbial sulfate reduction, probably via the oxidation of produced S2−98 and/or electron shuttling of produced sulfur intermediates94 (Fig. 7). In the study area, bacteria capable of disproportioning intermediate sulphur species was also detected in groundwater of the transition area and the flat plain with relative abundance of 1.6%, including Desulfocapsa99, Desulfurivibrio100, and an unclassified member of Thermodesulfovibrionia101 family. Moreover, sulfur and oxygen isotope analyses confirmed that both sulfide and sulfur intermediates (i.e., elemental sulfur) produced from microbial sulfate reduction reduced the Fe(III) oxides102. Hence, the positive correlation between S2- and Δ18OPO4 values (r = 0.40, p < 0.05; Fig. 5d) indicated that more PO4 with an equilibrium δ18OPO4 value was released into groundwater via the reductive dissolution of Fe(III) oxides, as a result of increasing microbial sulfate reduction under anoxic conditions (Fig. 7). Notably, both the concentrations of S2- and PO4 decreased from the transition area to the flat plain (Fig. 3f, h). This was possibly due to the precipitation of ferrous sulfide minerals, such as pyrite103, and subsequent PO4 adsorption104. In the flat plain, pyrite was oversaturated in groundwater37 and observed in aquifer sediments103. Since previous experiments demonstrated a negligible isotopic fractionation between dissolved and adsorbed PO4 during the abiotic process of PO4 adsorption on Fe-minerals at long reaction time ( > 100 h)64, we do not expect a shift in δ18OPO4 signals by PO4 adsorption on ferrous sulfide minerals in this study.

In addition, the increase in Δ18OPO4 value towards zero from the alluvial to the flat plain was related to the increasing groundwater residence time along the flow path during which dissolved PO4 in groundwater was more intensively cycled by microorganisms. In the study area, the average groundwater flow rates gradually decrease from 20 m yr−1 in the alluvial fan, over 10 m yr−1 in the transition area, to 5 m yr−1 in the flat plain105. Accordingly, we estimated that the mean residence times of groundwater in the alluvial fan area (around 2.5 km-length), in the transition area (around 1.5 km-length), and in the flat plain area (around 2 km-length) were 125 years, 275 years, and 675 years, respectively. During these long groundwater residence times, the abiotic oxygen isotopic exchange between dissolved PO4 and ambient groundwater should be discussed. Previous studies predicted that the half-live for a complete abiotic equilibrium exchange between PO4 and water would be 88,000 yrs at 10 °C44,106, which was close to groundwater temperature in the study area (11.3 °C on average; see Supplementary Table S1). The predicted time was much longer than the groundwater residence time in the study area. Therefore, along the flow path, the changes of Δ18OPO4 values more likely resulted from microbially-mediated isotopic exchange. The longer residence time of groundwater allowed the dissolved PO4 to be more microbially cycled and consequently the δ18OPO4 values of groundwater gradually increased towards isotopic equilibrium. This was in line with previous studies, demonstrating that the long residence time in the lower catchments of the Yasu River (Japan) and the Redon River (France) fostered biological PO4 recycling and the equilibration of δ18O values of dissolved PO439,40.

Methods

Study area

The study area is located in the northwest of the Hetao Basin in the Inner Mongolia Autonomous Region, PR China (Fig. 1), between the Langshan Mountains and the Yellow River. The Quaternary sediments that fill the basin have both alluvial and lacustrine sources, which are mainly derived from the metamorphic complex of the Langshan Mountains and partly from fluvial deposits of the Yellow River61,62. There are two clay layers that are enriched in organic matter at the depths of approximately 5 m and 40 m below land surface (bls), overlying two sandy and silty aquifers and preventing vertical filtration61. Groundwater flows from the Langshan Mountains in the north-west in a south-eastern direction towards the centre of the Hetao Basin. Here, hydraulic conductivity decreases from about 9.0 m d−1 near the mountains to 3.8 m d−1 in the centre of the floodplain79.

Along the groundwater flow, there is a typical hydrogeochemical zonation in the aquifers. In the alluvial fan near the mountains, the aquifers consist of pluvial-alluvial sediments, whereas the aquifer sediments in the flat plain are mainly composed of inland lacustrine sediments61. A transition area was identified between the alluvial fan and the flat plain, which is characterized by mixed-source sediment deposits62. In the aquifer sediments of the study area, high PO4 contents were mostly found in clay/silty clay layers, whilst PO4 contents in sand/silt layers were comparatively low. Ca/Fe(II)-minerals represented the dominant sedimentary PO4 pool of total extractable PO4, followed by exchangeable PO4, PO4 in Fe(III) oxides, strongly adsorbed PO4 to metal-(oxyhydr)oxides and residual PO437. Moreover, along the groundwater flow path, Porg contents were obviously higher in the sediments of the flat plain than in the sediments of the alluvial fan37. Based on our previous investigation at the same site, the average groundwater flow rates gradually decrease from 20 m yr−1 in the alluvial fan, over 10 m yr−1 in the transition area, to 5 m yr−1 in the flat plain105. From the alluvial fan to the flat plain, distinct redox reactions occur with the sequence of O2-reduction, NO3-reduction, Mn(VI)-reduction, As(V)-reduction, Fe(III)-reduction and SO42-reduction59. Representative microbial communities that catalyze those redox reactions were also detected, including Gallionella (Fe(II)-oxidizer), Pseudomonas (denitrifier), Rhodoferax (Fe(III)-reducer), Bacillus (As(V)-reducer) and Desulfatiferula (sulphate-reducer)60.

Groundwater sampling and field analysis

Along the groundwater flow path, 45 groundwater samples, including 13 samples from the alluvial fan, 16 samples from the transition area, and 16 samples from the flat plain, were collected from the irrigation wells with depths of 10–100 m bls in July 2019 (Fig. 1). For data evaluation, groundwater samples from shallow and deep aquifers were pooled for the respective zones. Based on the previous investigation, the spatial variations of physiochemical groundwater parameters, including PO4 concentrations37, showed no significant difference (p > 0.05) between the two aquifers along the groundwater flow path. Hence, we did not differentiate the shallow and deep aquifers in this study.

Before sampling, wells were pumped for more than 20 min until sensitive parameters, including groundwater temperature (T), electrical conductivity (EC), pH, oxidation reduction potential (ORP, which was converted into Eh after corrected to the standard hydrogen electrode provided in the manual of HANNA, HI 9828), and dissolved oxygen (DO) of groundwater were stable. All samples used for laboratory analysis were filtered through 0.45 μm membrane filters (cellulose nitrate filters, Jinteng), with the exception of samples used for dissolved organic carbon (DOC), which were filtered through 0.45 μm syringe filters (PTFE filters, Jinteng). All bottles and flasks for sample storage were acid-washed and preconditioned before use. Details on groundwater sampling and field analysis please see Supplementary Note 1.

For groundwater δ18OPO4 analysis, a volume of 1–3 L of fresh groundwater was collected in polyethylene bottles to obtain sufficient PO4 for the isotopic analysis (around 1 mg in total), depending on the PO4 concentrations measured in the field (ascorbic acid; PhosVer3; Hach, DR2800). To isolate PO4 from groundwater, we added 1–3 g of FeSO4 (Sigma-Aldrich, ACS reagent, ≥ 99.0%) into the groundwater to advance the coprecipitation of PO4 by Fe(III)-oxides after 24 h of aeration107. The precipitates were filtered and air-dried, and then stored in sterile polythene bottles for further laboratory analysis. Phosphate concentrations in filtrates were analyzed after the filtration to make sure that PO4 in groundwater was completely removed by the co-precipitation with Fe(III)-oxides. Results showed that PO4 concentrations in filtrates taken after co-precipitation and filtration remained below the detection limit (5 μg L−1 for UV-VIS spectrophotometer [Varian Cary 100] based on the molybdenum blue method108), indicating a successful separation of PO4 from groundwater.

Sediment sampling and PO4 extraction

Two boreholes were both drilled with depths between 10 and 80 m below the land surface (bls). Specifically, borehole K2 was located in the alluvial fan area near the mountains, while borehole K1 was located in the flat plain. Four sediment samples were collected from each borehole with increasing depths (see Supplementary Table S2). After taken out, sediment samples were immediately sealed in foil and in sterile bags flushed with ultrapure nitrogen gas. Before the laboratory PO4 extraction and δ18OPO4 analysis, samples were preserved in a refrigerator at −80 °C.

To differentiate PO4 originating from microbial cycling and from dissolution of primary minerals, we targeted the PO4 in apatite minerals109,110. In order to extract HCl-PO4, we first used NaOH solution to remove labile P (including exchangeable P, adsorbed P onto Fe/Mn-(oxyhydr)oxides and structurally bound P in Al-(oxyhydr)oxides)7. Then, HCl-PO4 in sediments was extracted by 1 M HCl solution (Sigma-Aldrich, p.a. ≥ 37 %)37. According to previous studies, PO4 extracted by HCl solution could mainly contain three different sedimentary P pools, including Fe(III) oxides-bound PO4, detrital apatite and authigenic apatite44,111,112. Prior to HCl-PO4 extraction, we had checked the presence of Fe(III) oxides-bound PO4 by the extraction of citrate-dithionite-bicarbonate (CDB) solution37. The results showed that Fe(III) oxides-bound PO4 (average content: 9.62 ± 2.65 mg kg−1; n = 8) only accounted for a negligible proportion (about 4%) if it was co-extracted with HCl-PO4 (average content: 242 ± 24.7 mg kg−1; n = 8). Moreover, adsorbed PO4 onto Fe(III)-(oxyhydr)oxides was removed by NaOH solution7. Hence, HCl-PO4 in this study mainly represented sedimentary P pool of apatite minerals, including detrital apatite and authigenic apatite. Details on extraction of HCl-PO4 from sediments please see Supplementary Note 2.

Laboratory analysis

Stable oxygen isotope ratio analyses of PO418OPO4) were carried out using groundwater and sediment samples from the Hetao Basin. Prior to the δ18OPO4 analyses of groundwater samples and extracted HCl-PO4 from sediments, a multistep purification procedure was modified based on the protocols developed by Tamburini et al.110 and Neidhardt et al.15. In brief, this procedure consists of multiple dissolution, precipitation, and resin sorption steps to remove all other possible O-bearing compounds and constituents that would interfere with the isotopic analysis (see Supplementary Fig. S1). According to our previous investigation in the Hetao Basin, the average concentrations of dissolved PO4 and dissolved Porg in groundwater were 82.4 and 6.2 μg L−1, respectively, which suggested that dissolved P in groundwater mainly presents in the form of PO4. Moreover, Porg components in the sediments mostly preserves in the exchangeable P pool37, which was removed by NaOH solution prior to HCl extraction. Hence, the isotopic effect on δ18OPO4 analyses of groundwater and HCl-extractable PO4 by means of Porg hydrolysis is considered negligible in this study. Details on the laboratory analyses of groundwater major cations, anions, dissolved P species, DOC, and oxygen and hydrogen isotopes (δ18OH2O and δ2HH2O) please see Supplementary Note 3.

For PO4 purification from groundwater, the Fe(III)-(hydr)oxide precipitates were dissolved in 50 mL 16.5% HNO3 to release co-precipitated PO4 and the filters were rinsed with Millipore water to completely remove PO4. The extracted HCl-PO4 from sediments were preserved in 50 mL 1 M HCl solutions. Then, the protocol basically followed the method from reference 110, which consisted of two precipitation steps (see Supplementary Fig. S1). First, PO4 was precipitated in the form of ammonium-phospho-molybdate (APM, (NH4)3PMo12O40) by adding 25 mL 35% NH4NO3 solution (Roth) and 40 mL 10% ammonium molybdate solution ((NH4)6Mo7O24·4H2O) (Roth, 99%) and placing in a warm bath at 20°C. After filtration (0.22 μm PES membrane; Millipore Express, Merck) and rinsing with Millipore water, the APM precipitate was dissolved in 25 mL ammonium-citrate solution. The next precipitation step was then induced by adding 25 mL acidic magnesia solution and 6 mL NH4OH/H2O (1:1 v/v) solution in order to keep the pH between 8 and 9. Magnesium-ammonium-phosphate (MAP, NH4MgPO4·6H2O) crystals form while stirring overnight. After filtration (0.22 μm PES membrane; Millipore Express, Merck) and rinsing, MAP crystals were dissolved by adding 20 mL 0.5 M HNO3. After dissolution, possible interfering cations were removed by shaking the samples together with 6 mL cation exchange resin slurry (BioRad, AG® 50W-X8, H+ form) overnight. Resin and sample solution were then separated by filtration (0.22 μm PES membrane; Millipore Express, Merck). Afterwards, 4 mL 2 M AgNO3 was added to the sample solution to remove chloride as AgCl precipitate to assure the purity of the final Ag3PO4 precipitate. Sample solutions were filtered (0.2 μm GTTP membrane; Isopore, Merck) again and then transferred to new 50 mL tubes. Finally, the PO4 was precipitated as Ag3PO4 following the addition of 5 mL silver ammine solution and the placement of the open tubes into an oven at 50 °C for 48 h. The δ18OPO4 of groundwater and HCl-PO4 from sediments was then analysed in form of pure Ag3PO4 by TC/EA-IRMS (PYRO cube and IsoPrime100, Elementar Analysensysteme, Germany) in a continuous flow mode. The isotope analysis was carried out in the laboratory of Soil Science and Geoecology at the University of Tübingen. All δ18OPO4 values were reported in the conventional per mil notation (‰) relative to the Vienna Standard Mean Ocean Water (VSMOW) reference standard.

For internal calibration and quality control of groundwater δ18OPO4 measurements, three certified international reference standards (USGS80, δ18OPO4 = +13.1‰; USGS81, δ18OPO4 = +35.4‰; IAEA602, δ18O = +71.3‰) were included in each TC/EA IRMS run. Within the sample weight range of 0.1–0.7 mg, the standard error (SE) of the repeatedly measured δ18OPO4 values of the three reference standards were ± 0.07‰ (USGS80; n = 20), ± 0.26‰ (USGS81; n = 18), and ± 0.97‰ (IAEA602; n = 8), respectively (see Supplementary Table S3). The differences between the average measured δ18OPO4 values and certified values of the three reference standards were +0.2‰, −0.8‰, and +0.3‰, respectively (see Supplementary Table S3). Notably, the measured δ18OPO4 values of USGS80, which were the closest to those of the groundwater δ18OPO4 values, showed the smallest SE and the smallest difference between average values and certified values (see Supplementary Fig. S2a), reflecting a good stability and accuracy of the isotope analyses. All samples were analyzed as triplicates and within the same weight range as to the reference standards, except when the Ag3PO4 mass was too low. The oxygen yield was checked by comparing the weights of Ag3PO4 introduced into TC/EA-IRMS and the obtained CO peak area50. The purity of the precipitated Ag3PO4 was ensured by the close match of two regression lines representing the standards and the samples (see Supplementary Fig. S2b). Moreover, possible effects of arsenic contamination on the δ18OPO4 values were considered15. Results showed that arsenic was mostly removed from the samples through the procedure, with the effects on δ18OPO4 values being less than 1‰ (see Supplementary Table S4). The δ18OPO4 measurements of HCl-PO4 extracted from sediments were also calibrated by two internal standards (USGS80, δ18OPO4 = +13.1‰; USGS81, δ18OPO4 = +35.4‰) and the quality was controlled by an additional external standard (IAEA601, δ18O = +23.1‰). Three certified standards were included in each TC/EA IRMS run and the results also showed a good accuracy and reproducibility for all reference materials (see Supplementary Table S5).

Pure solid KH2PO4 salt (Roth, ≥ 99.0%) was dissolved and included as an internal quality control in each sample batch. Mass balance results showed that after the purification, an average of 89.3 ± 0.3% (n = 4) of the initial PO4 was converted into Ag3PO4 precipitate. The δ18OPO4 value of the processed internal standard shifted approximately 0.35‰ VSMOW (average: +13.98 ± 0.12‰, n = 12) when compared to that of the pure initial KH2PO4 salt (average: +13.63 ± 0.14‰, n = 3). Thus, the internal standard demonstrated that PO4 loss and oxygen fractionation during the purification were insignificant (see Supplementary Tables S6 and S7). Further details regarding the quality of the δ18OPO4 analyses of groundwater and sedimentary HCl-PO4 please see Supplementary Fig. S2 and Tables S3S7.

16S rRNA gene analysis

At the study site, 12 groundwater samples used for the δ18OPO4 analyses in this study have been previously analyzed for 16S rRNA gene amplicon sequencing60,113,114. They include 4 samples in the alluvial fan, 4 samples in the transition area and 4 samples in the flat plain. For all those 16S rRNA gene analyses at the study site, DNA was extracted from each sample using the same FastDNA SPIN Kit for Soil (MP Biotechnology, Solon, OH, USA) and the 16S rRNA gene was amplified from the genomic DNA using the same primers 338 F (ACTCCTACGGGAGGCAGCAG) and 806 R (GGACTACHVGGGTWTCTAAT) targeting the V3-V4 region60,113,114. In this study, raw FASTQ files containing results of 16S rRNA gene sequencing were collected from those previous studies. Then, the sequences were imported into QIIME2 v2023.5115 and processed with DADA2116, including trimming, quality-filtering and demultiplexing. The remaining sequences were clustered into operational taxonomic units (OTUs) at the 97% sequence identity level. Finally, the taxonomic assignment of each OTUs was conducted using Naïve Bayes classifier algorithm trained on data from SILVA v.138.1. The clean data for each sample were uploaded to the NCBI database under the accession number PRJNA1082075. More details regarding the DNA extraction and high-throughput sequencing please refer to Xiu et al.60, Ke et al.113 and Wang114.

Data processing

Microbial PO4 turnover that is catalyzed by intracellular enzymes (i.e., pyrophosphatases) leads to an exchange of O atoms between PO4 and ambient water (groundwater in our case) associated with a temperature-dependent equilibrium fractionation65. The resulting theoretical equilibrium value of δ18OPO418OPO4 equ) after microbial PO4 turnover can be calculated based on known temperature (T) and stable oxygen isotope value of ambient water (δ18OH2O)107. In the study area, despite the relatively long residence times of groundwater ( ~ 100 yr - ~1200 yr105) (see Study area in Methods), the relatively fast intracellular isotope fractionation can reach equilibrium within several days50, which is why we assume that measured groundwater δ18OH2O and temperatures were representative for the calculation of δ18OPO4 equ. Recently, the equation was revised by ref. 50, and we rearranged the δ18OPO4 equ calculation as Eq. (1).

$${{{\rm{\delta}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{PO}}}}4{{{\rm{equ}}}}}=\left({{{\rm{\delta}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{H}}}}2{{{\rm{O}}}}}+1000\right)* {{{{\rm{e}}}}}^{\left(14.43* 1000/{{{\rm{T}}}}[{{{\rm{K}}}}]-26.54\right)/1000}-1000$$
(1)

An isotopic shift of measured values of δ18OPO4 from the theoretical equilibrium value of δ18OPO4 equ18OPO4) is defined as Eq. (2)49.

$$\Delta ^{18}{{{{{\rm{O}}}}}_{{{{\rm{PO}}}}4}}={{{\rm{\delta}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{PO}}}}4}-{{{{\rm{\delta }}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{PO}}}}4{{{\rm{equ}}}}}$$
(2)

Based on the observed δ18OPO4 and calculated δ18OPO4 equ in groundwater, a simple two-end member mixing model was applied in this study to calculate the contribution of microbial cycling to dissolved PO4 in the aquifers, relative to other abiotic processes that release PO4 with specific O isotope values.

Extracellular Porg hydrolysis is catalyzed by enzymes and occurs in two steps, including (1) hydrolytic cleavage of phosphodiester (P-diester) bonds that results in the generation of phosphomonoesters (P-monoesters) and (2) conversion of P-monoesters to an organic group plus PO447. The first step is catalyzed by phosphodiesterase (diesterase) and the second step is catalyzed by phosphomonoesterase (monoesterase)47. For each step, PO4 inherits three atoms from the Porg molecule and incorporates one O atom from ambient water70,71. Hence, PO4 released from extracellular Porg hydrolysis contains signatures of both enzymatic steps and the combined isotopic effect of diesterase + monoesterase is expressed as Eq. (3):

$${{{\rm{\delta}}}} ^{18}{{{{{\rm{O}}}}}_{{{{\rm{PO}}}}4}}=0.5\times {{{\rm{\delta}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{org}}}}}+0.5\times \left({{{\rm{\delta}}}}^{18}{{{{\rm{O}}}}}_{{{{\rm{water}}}}}+{{{\rm{F}}}}\right)$$
(3)

where δ18Oorg stands for the original stable oxygen isotope value of organophosphorus compounds and F is the combined fractionation factor associated with breaking diester and monoester bonds (F = (F1 + F2)/2)47. In natural waters, PO4 is expected to have a stronger RNA degradation signature than DNA degradation signature47. For RNA, the fractionation factor between O incorporated into PO4 and ambient water O is −5( ± 6)‰ for the diesterase + alkaline phosphatase (APase) pathway (F1 = +20‰; F2 = −30‰)70,71.

Along the groundwater flow path, distributions of δ18OPO4 values and redox components in the three zones were analyzed by SPSS Statistics (Vers. 25) using group comparisons. The data were first checked for normal distribution using the Shapiro-Wilk test. Based on the outcome, group means of variables in different zones were compared by a one-way ANOVA, or a non-parametric Kruskal-Wallis one-way ANOVA on ranks and by a pair-wise post-hoc test (LSD) with Bonferroni correction. The group means of Δ18OPO4 in different zones were compared with zero by a one-sample t-test. Moreover, two-tailed Spearman correlation tests and ordinary least squares regression analysis were applied to examine the relationships between the measured variables. In addition, we evaluated the residuals of selected regressions. The confidence interval for all tests was set to 95% and p-values ≤ 0.05 were considered statistically significant.

Reporting summary

Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.