Skip to main content

Thank you for visiting You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

Carbon fixation rates in groundwater similar to those in oligotrophic marine systems


The terrestrial subsurface contains nearly all of Earth’s freshwater reserves and harbours the majority of our planet’s total prokaryotic biomass. Although genetic surveys suggest these organisms rely on in situ carbon fixation, rather than the photosynthetically derived organic carbon transported from surface environments, direct measurements of carbon fixation in the subsurface are absent. Using an ultra-low level 14C-labelling technique, we estimate in situ carbon fixation rates in a carbonate aquifer. We find these rates are similar to those measured in oligotrophic marine surface waters and up to six-fold greater than those observed in the lower euphotic zone. Our empirical carbon fixation rates agree with nitrification rate data. Metagenomic analyses reveal abundant putative chemolithoautotrophic members of an uncharacterized order of Nitrospiria that may be behind the carbon fixation. On the basis of our determined carbon fixation rates, we conservatively extrapolate global primary production in carbonate groundwaters (10% of global reserves) to be 0.11 Pg carbon per year. These rates fall within the range found for oligotrophic marine surface waters, indicating a substantial contribution of in situ primary production to subsurface ecosystem processes. We further suggest that, just as phototrophs are for marine biogeochemical cycling, such subsurface carbon fixation is potentially foundational to subsurface trophic webs.


The continental subsurface is the planet’s largest carbon reservoir1, housing up to 19% of its total biomass2,3 and 99% of its freshwater4. Despite accounting for only 6% of total stores, modern groundwater (the fraction accrued in aquifers over the past 50 years) is the single most significant source of potable water. Carbonate karst aquifers alone are thought to supply people with nearly 10% of their drinking water5. Unfortunately, modern groundwater is also the most vulnerable to anthropogenic and climatic impacts4. While subsurface ecosystems have long fascinated ecologists6, and more recently microbiologists7, accessibility, enormous spatial heterogeneity and challenges in interpreting process rate measurements have obscured a meaningful understanding of their contributions to global biogeochemical cycles8.

The widespread recognition that Earth’s biosphere extends deep into the subsurface occurred only recently9. Historically, carbon supply in such environments was thought to be limited to the trickling of surface-produced organic matter into the shallow subsurface10 or what was stored within sedimentary rocks11. By stark contrast, a wealth of compelling genetic evidence suggests that in situ carbon fixation is critical for sustaining highly diverse microbial metabolic networks in groundwater, in both the shallow and the deep subsurface12,13,14,15,16,17,18,19. Despite the implications of gene-based surveys, the empirically derived activity measurements required to corroborate such inferences, constrain biogeochemical fluxes, understand system dynamics and integrate processes into regional and global models have yet to be reported. Here we report our use of a novel radiocarbon method to derive empirical carbon fixation rates and place them in the context of global groundwater.

Groundwater CO2 fixation rates resemble marine waters

Primary productivity in the shallow subsurface (groundwater wells 5–90 m deep), experimental carbon fixation rates varied from 0.043 ± 0.01 to 0.23 ± 0.10 μgC l−1 d−1 (mean ± 1 standard deviation (s.d.); Fig. 1a, Extended Data Table 1 and Supplementary Information). The ultra-low level 14C-labelling approach developed in this investigation exploits the high sensitivity of accelerator mass spectrometry, thereby minimizing impacts to groundwater hydrochemical equilibria and affording shorter incubation times. This method is particularly useful within a carbonate geological setting, where high dissolved inorganic carbon (DIC) backgrounds and a scarcity of microbes warrant greater sensitivity than is achievable via scintillation-based 14C-labelling approaches. Rates resulting from our new labelling technique probably approximate net primary productivity rather than gross productivity, as has been reported for marine systems20,21, and we further expect them to be conservative estimates for carbon fixation as they consider only contributions from the planktonic portion of the community (Supplementary Information).

Fig. 1: Carbon fixation rates and metagenomic potentials across the aquifer.
figure 1

a, Rates of carbon fixation. Outer error bars depict one standard deviation; inner grey bars delineate standard error of the mean. Rates for well H51 are derived from non-labelled controls (Supplementary Information). Letters denote the results of ANOVA and post hoc Tukey tests (excluding H51). b, Relative importance of the predicted carbon fixation pathways. c, Electron donor sources in each well. Values are averages from the triplicate 0.2 μm filtered fraction metagenome samples. d, Mean dissolved oxygen (DO) concentrations in groundwater collected in summer months (May–September 2010–2018); error bars depict standard deviation. e, Redox-potential measurements from identical time points; error bars depict standard deviation. While groundwater from H52 and H43 exhibited anoxic or hypoxic conditions, the positive redox potentials were due to oligotrophic conditions and mean nitrate concentrations of ~4 mg l–1.

We compared these carbon fixation rates that were measured in groundwater of varying biogeochemical characteristics22 with the only other subsurface 14CO2 assimilation measurements reported: those of a deep (830–1,078 m) groundwater borehole from crystalline bedrock in Sweden23. To do so, we converted the published rates of isotopic incorporation to carbon equivalents, revealing the lower but overlapping range of 0.0095 to 0.0560 μgC l−1 d−1.

To better understand the relevance of the rates measured, we compared them with those of well-documented oligotrophic marine surface waters. Unlike our samples, the carbon fixed in these waters was sourced almost entirely by bacterial photoautotrophs24,25. When compared directly with a comprehensive dataset compiled by ref. 26, our rates overlapped with those of global marine waters at depths to 140 m, equating to roughly 10% of the reported global median for 0–20 m depths (2.65 μg l−1 d−1, interquartile range (IQR) = 1.74, 6.02) and 20% of the median for 20–140 m depths (1.2 μg l−1 d−1, IQR = 0.6, 1.7). Comparisons with the extensively studied Sargasso Sea in the Bermuda Atlantic Timeseries Study (BATS)27,28 and the Hawaiian Oceanographic Timeseries (HOTs)29 datasets yielded similar findings (Fig. 2). Our rate measurements ranged between 3% and 23% of the median reported net primary productivity in the upper euphotic zones (down to ~120 m) and between 20% and 600% of the median of the lower euphotic region (100–120 m).

Fig. 2: Comparison of carbon fixation rates within groundwater and the marine euphotic zone.
figure 2

Violin plots depicting the distribution of carbon fixation rates measured in oligotrophic marine surface waters and groundwater. HOTs, Hawaiian Oceanographic Timeseries (1999, cruises 101–110); BATs, data from 1988 to 2016 for the Bermuda Atlantic Timeseries; Liang, a collated dataset compiled by Liang et al.26; GW, the range of groundwater samples shown in Fig. 1.

We also considered contributions to existing particulate organic carbon (POC) stocks and new carbon inputs per microbial cell count. After normalizing for estimated total bacterial cell numbers, groundwater yielded 0.3–10.8 fg fixed carbon per bacterial cell per day (Extended Data Table 2), which matched estimates of 0.25–12.1 fgC per bacterial cell per day across the marine photic zone (5–150 m). However, groundwater received new daily carbon inputs of only 0.47% ± 0.22% (mean ± s.d.; Extended Data Table 2) of its existing POC, much lower than the marine system’s 2.6% ± 2.9% gain in the lower euphotic zone and 22% ± 18% at the surface30,31. This disparity might stem from the larger recalcitrant fraction of POC in groundwater compared with oligotrophic oceans, which is supported by deviations in 14C and 13C signatures of total groundwater POC concentrations compared with lipid signatures of resident microbes32.

An ecosystem dominated by chemolithoautotrophs

To identify dominant microbial primary producers, a dereplicated and quality-controlled set of 1,224 metagenome-assembled genomes (MAGs) were generated from groundwater samples. Of these, 102 putative chemolithoautotrophs exhibiting at least 50% completion scores for carbon fixation pathways were identified (Fig. 3). Almost exclusively bacterial (101), these MAGs represented 17 distinct phyla, 21 classes and 35 families (Fig. 3 and Supplementary Dataset 1). In some samples, up to 12% of metagenomic reads from a sample could be recruited to these potentially chemolithoautotrophic MAGs (Supplementary Information and Extended Data Figs. 1 and 2). A single archaeal MAG of the family Nitrosopumilaceae encoded gene products for the 4-hydroxybutyrate/3-hydroxypropionate pathway and was not relatively abundant (<5× maximum normalized coverage).

Fig. 3: Phylogeny, relative abundances and transcriptional activities of putative chemolithoautotrophic MAGs.
figure 3

Approximately maximum-likelihood phylogenetic tree based on concatenated single-copy protein alignments for all bacterial MAGs considered. Branches are coloured according to the predicted carbon fixation pathway, and the matching leaf is indicated by a point. Bar charts present average normalized metagenomic coverages within each well from triplicate 0.2 μm filtered fractions. Pie charts show the coverage of mRNA transcripts recruited, normalized to gene number and library size. Asterisks denote MAGs discussed in greater detail in Supplementary Information. The tree is rooted using Patescibacteria (CPR) as an outgroup, indicated by the collapsed grey leaf in the upper left. Node numbers: (1) c_Nitrospiria, (2) c_Thermodesulfovibrionia, (3) o_Nitrospirales, (4) f_Nitrospiraceae, (5) g_Nitrospira_D, (6) g_Nitrospira, (7) p_Nitrospinota, (8) c_Gammaproteobacteria, (9) o_Acidiferrobacterales, (10) f_Sulfurifustaceae, (11) g_SM1-46, (12) f_UBA6901, (13) o_Burkholderiales, (14) f_Nitrosomonadaceae, (15) g_Nitrosomonas, (16) f_SG8-41, (17) f_SG8-39, (18) c_Brocadiae, (19) o_Brocadiales, (20) f_Brocadiaceae and (21) f_Scalinduaceae.

Three chemolithoautotrophic pathways were detected (Fig. 1b): the Calvin–Benson–Bassham (CBB), Wood–Ljungdahl (WL) and reverse tricarboxylic acid (rTCA) cycles were present in 37, 50 and 15 MAGs, respectively. The summed and normalized relative coverages of MAGs equipped with these metabolic pathways aligned with the carbon fixation rates measured in wells H52, H32 and H14 while contrasting with rate data from wells H41 and H43 (Fig. 1, Extended Data Fig. 2 and Supplementary Information). The greatest relative abundances of predicted chemolithoautotrophs were detected in oxic well H41 and anoxic well H52. Anoxic groundwater was dominated by putative sulfur-oxidizing (53% of summed and normalized coverages of all chemolithoautotrophic MAGs) and putative anaerobic ammonium-oxidizing (anammox (10%)) autotrophic microbes, while oxic groundwaters harboured greater abundances of potential nitrifiers (76%; Fig. 1c and Supplementary Information).

Uncharacterized microbes influence CO2 fixation potential

The most abundant putative chemolithoautotrophic populations represented by MAGs generated from anoxic groundwater were of poorly studied and/or uncharacterized microbial lineages. Those most abundant in oxic groundwaters, however, were phylogenetically and metabolically similar to well-characterized microbes (Supplementary Information). In both cases, metabolic reconstructions suggested that dominant subpopulations could access a diverse suite of (in)organic electron acceptors and donors. We mapped previously generated transcripts to these MAGs to confirm the active expression of gene products involved in energy acquisition and carbon fixation. As opposed to the broad distributions posited by deoxyribonucleic acid- (DNA-) based abundances, transcript data revealed far more restrictive ranges in which specific gene products were favoured (Fig. 3). Given their metabolic versatility and the results of previous cultivation-based analyses33, these populations are expected to be mixotrophic (capable of supplementing carbon requirements with available organic matter). Overall, carbon fixation in anoxic groundwater was predicted to be fuelled by reduced sulfur, and there were three highly abundant, putative sulfur-oxidizing MAGs identified, each accounting for >2% of the total metagenomic reads in some samples (100–400× normalized coverages). Diverse reduced sulfur species fuelling these metabolisms are released from pyrite weathering of the karst rock34.

The most abundant MAG (H51-bin250-1) encountered in this study belongs to a deep-branching order, 9FT-COMBO-42-15, of class Nitrospira and is the first representative of class Nitrospiria thought to fix carbon via the WL pathway (Fig. 3, Fig. 4a and Supplementary Information). As there is precedence for autotrophic WL-utilizing bacteria within phylum Nitrospirota, and Candidatus Magnetobacterium was characterized with an equally flexible metabolism, these traits may be more widespread within the phylum than previously thought35. In addition, two MAGs with the potential to couple sulfur oxidation to carbon fixation via the CBB cycle were identified as members of the Sulfurifustaceae family of Proteobacteria (Supplementary Information and Fig. 4b; H32-bin014, H32-bin069). These MAGs recruited tenfold more transcripts than their Nitrospirota counterparts and were among the most transcriptionally active putative chemolithoautotrophic genomes detected (Figs. 3 and 4, Extended Data Fig. 4 and Supplementary Information). With its closest reference genomes Sulfuricaulis limicola and Ca. Muproteobacteria (RIFCSPHIGHO2_12_FULL_60_33), the taxonomic identity of this family is under debate. Per Genome Taxonomy Database (GTDB) classification nomenclature, Muproteobacteria belong to the Sulfurifustaceae family, and members of this family have been posited to oxidize sulfur in both aquatic and terrestrial environments36,37,38.

Fig. 4: Metabolic reconstructions of dominant putatively chemolithoautotrophic MAGs.
figure 4

Bar charts below each metabolic model summarize the average normalized coverage across each sample, scaled proportionally. Values within the bar chart indicate the sum normalized coverage for each MAG. Balloon plots depict normalized transcript coverages for genes affiliated with each pathway. If multiple copies were present, only the most active copy was plotted. The text information over each panel includes the predicted taxonomy (a, Nitrospiria; b, Sulfurifustaceae; c, Brocadiaceae), MAG identifier and estimated % completion/% redundancy. DNRA, dissimilatory nitrate reduction to ammonia.

Planctomycetota MAGs, predicted to couple anaerobic ammonium oxidation to carbon fixation via the WL pathway, exhibited mean transcriptional activities on par with their Sulfurifustaceae MAG counterparts (Figs. 3 and 4c, Extended Data Fig. 5 and Supplementary Information). The elevated transcriptional activity of gene products within the CBB and WL pathways suggests that these taxa play a disproportionately large role in chemolithoautotrophy relative to their DNA-based abundances. Surprisingly, all putative anammox MAGs detected were transcriptionally active in oxic groundwater (Figs. 3 and 4b; wells H41 and H51). Anammox reactions are typically inhibited in the presence of oxygen39, although microbes will still express critically important genes in low oxygen environments40,41.

N-based rate measurements validate carbon fixation rates

To evaluate the relationship between anammox and carbon fixation in anoxic groundwaters, we compared the rates of each in a well harbouring the greatest relative abundance of putative anaerobic ammonium-oxidizing bacteria (well H52). Here, anammox rates of 1.2 ± 0.5 nmol l−1 d−1 N2 were measured, similar to rates in another freshwater aquifer42. Empirical stoichiometric data demonstrate that 1.02 moles of N2 is produced via anaerobic ammonium oxidation for every 0.066 moles of CH2O0.5N0.15 reduced to biomass43. Assuming equivalent stoichiometry, the rate of carbon fixation via anammox in groundwater would be 0.93 ± 0.39 ngC l−1 d−1, more than 200 times lower than the 220 ngC l−1 d−1 measured. This result is corroborated by metagenomic data that suggest the high rate of carbon fixation in anoxic groundwater is more likely driven by reduced sulfur than by reduced nitrogen.

Metagenomic and metatranscriptomic data predicted that nearly all the organic carbon produced under oxic conditions in well H41 would be coupled to nitrification. To test this, we monitored the rate of aerobic ammonium oxidation in this well and recorded a mean production of 125.8 ± 5.9 nmol NO2– + NO3 l−1 d−1. Since the most abundant nitrifiers detected were most closely related to complete ammonium-oxidizing bacteria (Supplementary Information), we based our calculations on the 394 mg protein per mol of ammonia growth yields of Nitrospira inopinata, a comammox organism44. Assuming a cellular composition of C5H7O2N (ref. 45) and 55% protein content, we estimated a rate of 48.5 ± 1.9 ngC l−1 d−1, which was well within the range of error for our measured rate of 43 ± 13 ngC l−1 d−1 and confirms the importance of nitrification for carbon fixation at this site. Furthermore, the stoichiometry determined for oligotrophic marine rTCA nitrifiers46 of 0.0216 mol C/mol N matched our calculated ratio of 0.0276 ± 0.0084, indicating they are responsible for most of the fixed carbon.

Global estimates for groundwater primary productivity

There are an estimated 22.6 million km3 of groundwater on Earth4, 2.26 and 12.66 million km3 of which are housed in carbonate and crystalline aquifers, respectively. If we assume that our average rates accurately represent carbonate groundwater systems, then 0.108 ± 0.069 PgC (mean ± s.d.) is fixed every year in this global ecosystem (Extended Data Table 3). If the values reported from crystalline aquifers23 are representative of this environment, then another 0.15 ± 0.11 PgC would be fixed there annually. Collectively, the net primary productivity of ~66% of the planet’s groundwater reservoirs would total 0.26 PgC yr−1, approximately 0.5% that of marine systems and 0.25% of global NPP estimates47. As these projections exclude the missing contributions from groundwaters within siliciclastic and volcanic geologic settings and activities of attached microorganisms, global contributions to the carbon cycle are expected be many-fold higher.

We showed that conservative estimates of carbon fixation rates in a carbonate aquifer reached 10% of the median rates reported in oligotrophic marine surface waters and six-fold greater than those observed in the lower euphotic zone. Within oxic groundwaters, our carbon fixation method was independently validated by nitrification rate measurements. Normalizing carbon fixation rates by estimated bacterial numbers revealed equivalent carbon input (0.3–12 fgC per cell) for both marine and groundwater systems, despite the fact that daily inputs of new POC were 40 times greater in marine waters. This disparity makes sense since trophic webs are simpler in the subsurface, and the export of organic matter is constrained by long water residence times within the aquifer. Complementary metagenomic analyses revealed that groundwater carbon fixation is not dominated by a single functional guild but rather has contributions from diverse pathways and versatile microorganisms that are setting specific. As the majority of photosynthetically derived carbon in marine systems is labile (half-life <1 day), the findings of this study solicit new hypotheses regarding carbon cycling in the subsurface, particularly those positing newly synthesized carbon rather than surface-derived organic matter as the primary source of fuel for microbiota. Indeed, subsurface primary producers need to be considered as important to ecosystem processes as marine phototrophs are known to be in the surface ocean. Applying these rates of carbon fixation to ecosystem processes alters the way we think about these environments, challenges the importance of surface-derived organic matter fluxes on shallow subsurface functioning and establishes a framework broadly applicable across groundwater systems.


Site description

Groundwater samples were sampled from the Hainich Critical Zone Exploratory (NW Thuringia, Germany)22,48,49. This aquifer assemblage consists of a multistory fractured system composed of alternating layers of limestone and mudstone that developed along a hillslope of Upper Muschelkalk bedrock22. The primary aquifer, represented in this study by wells H41 and H51, is oxic and lies within the Trochitenkalk Formation (moTK). Primarily suboxic to anoxic, mudstone-dominated overhanging strata lie within the Meissner Formation (moM) and are represented here by wells H14 (moM—substory 1), H32 (moM—5, 6, 7), H43 (moM—8) and H52 (moM—3, 4). Geochemically, H32 and H41 coalesce into a single cluster while each of the other wells represents a distinct regime. Consistent with previous microbiological characterizations, however, each well studied represented a distinct community state50.

14C–DIC incorporation assay

This method, similar to a sensitive methane oxidation technique previously described51, is a modification of traditional 14C–CO2 primary productivity approaches52 predicated on the sensitivity offered by accelerator-based mass spectrometry. Groundwater was collected in July 2020 during sampling campaign PNK130, as described by ref. 19. After approximately three well volumes had been discharged and physicochemical parameters stabilized, groundwater was collected directly into nine pre-sterilized 2 l borosilicate bottles, from the bottom up. Bottles were then overfilled with greater than two volumes and sealed with gas-tight rubber stoppers. Triplicate samples from each well were then subjected to three treatments. A labelling treatment consisted of 6.77 × 10−7 mmol C–NaHCO3 that contained 200 Bq of activity (50 μCi; American Radiolabeled Chemicals) diluted to 9.38 Bq μl–1 with sterilized milliQ water, adjusted to pH 10 and verified using a scintillation counter. An advantage of this 14C technique is that the small amount of tracer added (representing 0.000006% of the total DIC) did not change the substrate concentration or influence conditions such as pH that could affect microbial populations. Kill controls were prepared in the same way, except 10 ml 50% ZnCl2 (w/v; final concentration 36.7 mM) was added to inhibit microbial activity. Unamended groundwater was also used as a control. All bottles were incubated in the dark at near in situ temperature for ~24 hours. Entire volumes were acidified to pH 4 with 3 M HCl, bubbled with N2 for one hour to remove DIC and then filtered through pre-baked (550 °C, eight hours) quartz fibres (47 mm, 0.3 um pore size, Macherey–Nagel QN-10) using pre-baked filter stands (EMD Millipore).

Filters were vacuum dried, sealed in quartz tubes with cupric oxide wire under vacuum and combusted at 900 °C for two hours. Evolved CO2 was purified cryogenically, measured as pressure in a known volume to determine C content and reduced to graphite for measurement by accelerator mass spectrometry at the WM Keck Carbon Cycle Accelerator Mass Spectrometry facility53. From the label incorporation and amount of carbon retained on the filters (Supplementary Data File 2), fixation rates were calculated using equation (1):

The technical variation was at most 3.6% (median = 0.78%) of the biological variation for the 14C measurements and was not considered in standard error of the mean calculations. Standard error of the mean was determined for both the 14C-based measurements (difference between two sets of triplicates, label and control, or label and kill controls) and POC measurements (all nine bottles from each well), separately. These errors were then propagated to yield the final error estimations. Analyses of variance and post hoc Tukey honestly significant difference (HSD) tests were conducted on resulting summary statistics (mean ± s.e.m.) using the following utility54. All 14C enrichment values were calculated using the differences between the 200 Bq-labelled samples and the 200 Bq-labelled kill controls. Rates calculated on the basis of no-label addition controls are presented in Extended Data Table 1. Data from global oligotrophic marine systems were included from Supplementary Data Sheet 126, the Bermuda Atlantic Timeseries years 1988–2016 via FTP27,55 and the Hawaiian Oceanographic Timeseries via FTP56. POC data from both sites were extracted from Dryad datasets generated by refs. 30,31. Bacterial cell number estimates for Hawaiian Oceanographic Timeseries were obtained from the FTP site57.

15N-isotope incubation experiments

Groundwater from wells H41 and H52 was collected in September 2018 and November 2018 to measure nitrification rates and anammox rates, respectively. Briefly, groundwater was collected into sterile glass bottles, from the bottom up, using a sterile tube. Bottles were then overfilled with three volume exchanges and sealed headspace free with silicone septa. Each sample was collected in triplicate alongside one control bottle per well. Samples were kept at 4 °C until they were processed (no more than 2 hours post-collection).

For nitrification measurements, 10 ml was removed from each sampling bottle (total volume 0.5 l) and replaced with N2 to analyse inorganic nitrogen and pH. Groundwater from control bottles was sterile filtered through a 0.2 µm filter (Supor, Pall Corporation). Sterile filtered 15N ammonium sulfate solution (98%, Cambridge Isotope Laboratories), serving as a substrate for ammonia-oxidizing prokaryotes, was then added to a final concentration of 50 µM. Samples were incubated at 15 °C in the dark sans agitation for five days. Ten-millilitre fractions were removed and replaced with N2 at the outset of the experiment and after 12, 24, 48, 70 and 120 hours via filtration through 0.2 µm filters; these fractions were stored at –20 °C for isotopic ratio mass spectrometry analyses. Additional 10 ml fractions were removed at intervals to monitor pH and inorganic nitrogen during the incubation.

For anammox rate measurements, sampling bottles (total volume 1 l) were flushed with N2 under sterile conditions for 30 minutes to remove all remnants of oxygen. Five-millilitre fractions were removed and replaced with N2 from each sample (and control) bottle to assess background 14NH4+ concentrations. Subsequently, samples were spiked with either (1) 50 µM 15NH4+ + 5 µM 14NO2 or (2) 5 µM 15NO2 as previously described58. Control bottles, serving as abiotic controls, were sterile filtered (0.2 µm filters; Supor, Pall Corporation) before flushing and the addition of nitrogen compounds. To facilitate destructive sampling at eight time points, groundwater (30 ml; in triplicate) was dispensed into sterile serum bottles leaving ~8 ml of headspace. Bottles were immediately sealed with butyl septa and crimp sealed and the headspace was purged with He. All bottles were then incubated in the dark at 15 °C sans agitation, and incubations were terminated after 0, 12, 24, 36, 48, 60, 72 and 96 hours by adding 300 µl 50% (v/w) aqueous zinc chloride solution.

Nitrification rates were determined on the basis of 15NO2 + 15NO3 production in incubations with 15NH4+. 15NO2- and 15NO3 were converted to N2 via cadmium reduction followed by a sulfamic acid addition59,60. The N2 produced (14N15N and 15N15N) was analysed on a gas chromatography isotope ratio mass spectrometer as previously described61. Rates were evaluated from the slope of the linear regression of 15N produced with time and corrected for the fraction of the NH4+ pool labelled in the initial substrate pool. The production of 15N-labelled N2 from anammox was analysed on the same isotopic ratio mass spectrometer as for nitrification rates and calculated as described62. Note, denitrification was not detected in any of the 15NO2 incubations. T tests were applied (P < 0.05) to assess whether rates were significantly different from zero (Extended Data Fig. 3).

DNA extraction and sample preparation

Samples used to generate metagenomic libraries were collected in January 2019 during sampling campaign PNK 110. For each sample replicate, approximately 50–100 l of groundwater was filtered sequentially through 0.2-µm- and 0.1-µm-pore-sized polytetrafluoroethylene (PTFE) filters (142 mm, Omnipore Membrane, Merck Millipore; Supplementary Data File 3). With the exception of H32 (which did not yield sufficient volumes), each well was sampled in triplicate. H32 was duplicated using a sample previously collected during campaign PNK108 (November 2018). Filters were frozen on dry ice and stored at –80 °C before extraction. DNA was extracted using a phenol-chloroform-based method, as previously described63, and resulting DNA extracts were purified using a Zymo DNA Clean & Concentrator kit. Metagenome libraries were generated with a NEBNext Ultra II FS DNA library preparation kit, in accordance with manufacturer’s protocols. DNA fragment sizes were estimated using an Agilent Bioanalyzer DNA 7500 instrument with High Sensitivity kits depending on DNA concentrations and recommendations of protocols (Supplementary Data File 3). Sequencing of the 32 samples was performed at the Core DNA Sequencing Facility of the Fritz Lipmann Institute using an Illumina NextSeq 500 system (2 × 150 bp). Resulting metagenomic library sizes ranged from 16.4 to 22.1 Gbp (mean = 19.6 Gbp; Supplementary Data File 3), and raw data were deposited into the European Nucleotide Archive under project PRJEB36523.

Metagenomic assembly and binning

Adaptors were trimmed and raw sequences subjected to quality control processing using BBduk v.38.5164. Assembly and binning were performed as previously described65. Briefly, all libraries were independently assembled into scaffolds using metaSPAdes v.3.1266, all of which were taxonomically classified per ref. 65. For individual assemblies, open reading frames (ORFs) were identified using Prodigal v.2.6.3 in meta mode67. To generate coverage profiles, all quality-assessed and quality-controlled (QAQC) sequences from each of the 32 metagenomic libraries were mapped back to each of the 32 scaffold databases using Bowtie2 v. in the sensitive mode68.

Scaffolds were binned using differential coverages and tetranucleotide frequencies with Maxbin269. In addition, ESOM and abawaca70 were used for both manual and automatic binning, based on tetranucleotide sequence signatures, using 3 kbp and 5 kbp or 5 kbp and 10 kbp as minimum scaffold sizes, respectively. DAS Tool71 was used with default parameters to reconcile resulting bin sets. Complete sets of bins from each of the samples were dereplicated using dRep v.2.4.072. All scaffolds, bin assignments, ORF predictions and taxonomic annotations were then imported into Anvi’o v.6.073. Each of the resulting 1,275 bins was manually curated in Anvi’o v.6, considering both coverage and sequence compositions. In the end, 1,224 bins passed the 30% completeness (median = 61%, IQR = (49%,73%)) and 10% redundancy (median = 0%, IQR = (0%,1.4%)) quality thresholds.

Characterizations of the MAGs

ORFs originating from all of the resulting MAGs were annotated using kofamscan74 with the ‘detail’ flag, and KO annotations were filtered using a custom script ( This utility preserves hits with scores of at least 80% of the kofamscan defined threshold, as well as those exhibiting a score >100 if there is no threshold. We elected to relax the default thresholds since all MAGs representing putatively chemolithoautotrophic microbes were verified manually, and we noticed that the best reciprocal blast hits with known reference sequences routinely scored below the kofamscan thresholds; that is, we favoured false positives over false negatives since we included a secondary verification step.

KEGGDecoder75 was used to assess the metabolic potential of five of the primary chemolithoautotrophic pathways: the CBB cycle, the WL pathway, the reverse citric acid cycle, the 4-hydroxybutyrate 3-hydroxypropionate pathway and the 3-hydroxypropionate bi-cycle. MAGs were examined in greater depth if a given pathway was >50% complete. The MAGs representing potential chemolithoautotrophs were re-annotated using the online BlastKoala server76 with essential steps verified through blast77 against the RefSeq database. A collection of HMM models was used to determine which form of Rubisco was detected, along with potential hydrogenases37. Using blastp77, dissimilatory bisulfite reductases (dsrAB) were compared with a database compiled by ref. 78 to predict whether the pathway operated in an oxidative or reductive manner. Blast was used to compare gene hits for narGH/nxrAB (nitrate reductase/nitrite oxidoreductase) with a custom database based on sequences presented within ref. 79.

All QAQC reads were remapped to a database consisting of only contigs of dereplicated MAGs. Normalized coverages for each of the MAGs were determined by scaling the resulting Anvi’o-determined coverages on the basis of the number of RNA polymerase B (rpoB) genes identified in the QAQC filtered reads. RpoB sequences were identified using ROCker with the precomputed model80. Scaling factors were calculated by dividing the maximum number of rpoB identified in the 32 metagenomic libraries by the number of rpoB detected in each sample. Reported values represent averages of the triplicates/replicates, unless stated otherwise. The taxonomy of each MAG was evaluated using the GTDB_TK tool kit81 in concert with the Genome Taxonomy Database (release 89)82,83 and its associated utilities67,84,85,86,87,88. Single-copy marker genes were identified and aligned with GTDB_TK for all bacterial MAGs, and a phylogenetic tree of the concatenated alignment was constructed using FastTree2 v.2.1.10 in accordance with the JTT + CAT evolutionary model. The resulting phylogenetic tree was then imported into iToL89 for visualization, and all MAGs were subjected to growth rate index analysis within each metagenomic library90.

Previously generated mRNA-enriched and post-processed metatranscriptomic libraries were procured from project PRJEB2878391. The groundwater source of these metatranscriptomes was collected in August and November 2015. QAQC filtered reads were mapped to MAGs using Bowtie2 v.2.3.5 in sensitive mode68, and the total number of rpoB transcripts from each metatranscriptomic library was determined, as described in the preceding for metagenomes. The transcriptomic coverages for each ORF from each MAG were determined using Anvi’o v.6 and normalized via scaling-factor calculations based on the total number of rpoB reads from the original metatranscriptome library (the coverage of each ORF from each MAG was normalized to a community-wide estimate of the transcriptional activity of a housekeeping gene in each sample). Means were determined considering all of the metatranscriptomes generated from a given well, including different sampling time points. While well H32 was sampled only once, mean values from all other wells account for three to four metatranscriptome coverages each. In addition, an average of the resulting normalized coverages for each MAG from each sample (sum of the MAG transcriptional coverage divided by the number of ORFs) was determined to estimate the relative transcriptional activity of the MAGs across the transect. Data were compiled and processed using R v.3.5.2 with Rstudio v.1.1.46392,93 and the tidyverse package94, and colour schemes were generated using the RColorBrewer utility95. All MAGs were deposited in project PRJEB36505’s data repository.

Data availability

The metagenomic raw data for this study was uploaded to the European Nucleotide Archive (ENA) under project PRJEB36505, while the individual sample assemblies and the metagenome-assembled-genomes were deposited into the ENA under project PRJEB36523. To improve ease of access, we also uploaded all MAGs and their associated data to the Open Science Framework (OSF) Project All raw and summarized accelerator mass spectrometry (AMS) data are available in Supplementary Data File 2, which was also uploaded to the same OSF project as the MAGs.


  1. Falkowski, P. et al. The global carbon cycle: a test of our knowledge of Earth as a system. Science 290, 291–296 (2000).

    Google Scholar 

  2. McMahon, S. & Parnell, J. Weighing the deep continental biosphere. FEMS Microbiol. Ecol. 87, 113–120 (2014).

    Google Scholar 

  3. Magnabosco, C. et al. The biomass and biodiversity of the continental subsurface. Nat. Geosci. 11, 707–717 (2018).

    Google Scholar 

  4. Gleeson, T., Befus, K. M., Jasechko, S., Luijendijk, E. & Cardenas, M. B. The global volume and distribution of modern groundwater. Nat. Geosci. 9, 161–167 (2016).

    Google Scholar 

  5. Stevanović, Z. Karst waters in potable water supply: a global scale overview. Environ. Earth Sci. 78, 662 (2019).

    Google Scholar 

  6. Poulson, T. L. & White, W. B. The cave environment. Science 165, 971–981 (1969).

    Google Scholar 

  7. Rusterholtz, K. J. & Mallory, L. M. Density, activity, and diversity of bacteria indigenous to a karstic aquifer. Microb. Ecol. 28, 79–99 (1994).

    Google Scholar 

  8. Smith, H. J. et al. Impact of hydrologic boundaries on microbial planktonic and biofilm communities in shallow terrestrial subsurface environments. FEMS Microbiol. Ecol. 94, fiy191 (2018).

    Google Scholar 

  9. Alexander, M. Introduction to Soil Microbiology (Wiley, 1977).

  10. Griebler, C. & Lueders, T. Microbial biodiversity in groundwater ecosystems. Freshw. Biol. 54, 649–677 (2009).

    Google Scholar 

  11. Krumholz, L. R., McKinley, J. P., Ulrich, G. A. & Suflita, J. M. Confined subsurface microbial communities in Cretaceous rock. Nature 386, 64–66 (1997).

    Google Scholar 

  12. Probst, A. J. et al. Differential depth distribution of microbial function and putative symbionts through sediment-hosted aquifers in the deep terrestrial subsurface. Nat. Microbiol. 3, 328–336 (2018).

    Google Scholar 

  13. Magnabosco, C. et al. A metagenomic window into carbon metabolism at 3 km depth in Precambrian continental crust. ISME J. 10, 730–741 (2016).

    Google Scholar 

  14. Stevens, T. O. & McKinley, J. P. Lithoautotrophic microbial ecosystems in deep basalt aquifers. Science 270, 450–455 (1995).

    Google Scholar 

  15. Tiago, I. & Veríssimo, A. Microbial and functional diversity of a subterrestrial high pH groundwater associated to serpentinization. Environ. Microbiol. 15, 1687–1706 (2013).

    Google Scholar 

  16. Mccollom, T. M. & Amend, J. P. A thermodynamic assessment of energy requirements for biomass synthesis by chemolithoautotrophic micro-organisms in oxic and anoxic environments. Geobiology 3, 135–144 (2005).

    Google Scholar 

  17. Momper, L., Jungbluth, S. P., Lee, M. D. & Amend, J. P. Energy and carbon metabolisms in a deep terrestrial subsurface fluid microbial community. ISME J. 11, 2319–2333 (2017).

    Google Scholar 

  18. Jewell, T. N. M., Karaoz, U., Brodie, E. L., Williams, K. H. & Beller, H. R. Metatranscriptomic evidence of pervasive and diverse chemolithoautotrophy relevant to C, S, N and Fe cycling in a shallow alluvial aquifer. ISME J. 10, 2106–2117 (2016).

    Google Scholar 

  19. Herrmann, M., Rusznyák, A. & Akob, D. M. Large fractions of CO2-fixing microorganisms in pristine limestone aquifers appear to be involved in the oxidation of reduced sulfur and nitrogen compounds. Appl. Environ. Microbiol. 81, 2384–2394 (2015).

  20. Peterson, B. J. Aquatic primary productivity and the 14C–CO2 method: a history of the productivity problem. Annu. Rev. Ecol. Syst. 11, 359–385 (1980).

    Google Scholar 

  21. Viviani, D. A., Karl, D. M. & Church, M. J. Variability in photosynthetic production of dissolved and particulate organic carbon in the North Pacific Subtropical Gyre. Front. Mar. Sci. 2, 73 (2015).

    Google Scholar 

  22. Kohlhepp, B. et al. Aquifer configuration and geostructural links control the groundwater quality in thin-bedded carbonate–siliciclastic alternations of the Hainich CZE, central Germany. Hydrol. Earth Syst. Sci. 21, 6091–6116 (2017).

    Google Scholar 

  23. Pedersen, K. & Ekendahl, S. Assimilation of CO2 and introduced organic compounds by bacterial communities in groundwater from southeastern Sweden deep crystalline bedrock. Microb. Ecol. 23, 1–14 (1992).

    Google Scholar 

  24. Partensky, F. & Garczarek, L. Prochlorococcus: advantages and limits of minimalism. Ann. Rev. Mar. Sci. 2, 305–331 (2010).

    Google Scholar 

  25. Karl, D. M., Hebel, D. V., Björkman, K. & Letelier, R. M. The role of dissolved organic matter release in the productivity of the oligotrophic North Pacific Ocean. Limnol. Oceanogr. 43, 1270–1286 (1998).

    Google Scholar 

  26. Liang, Y. et al. Estimating primary production of picophytoplankton using the carbon-based ocean productivity model: a preliminary study. Front. Microbiol. 8, 1926 (2017).

    Google Scholar 

  27. Steinberg, D. K. et al. Overview of the US JGOFS Bermuda Atlantic Time-series Study (BATS): a decade-scale look at ocean biology and biogeochemistry. Deep Sea Res. 2 48, 1405–1447 (2001).

    Google Scholar 

  28. Gundersen, K., Orcutt, K. M., Purdie, D. A., Michaels, A. F. & Knap, A. H. Particulate organic carbon mass distribution at the Bermuda Atlantic Time-series Study (BATS) site. Deep Sea Res. 2 48, 1697–1718 (2001).

    Google Scholar 

  29. Karl, D. M. & Lukas, R. The Hawaii Ocean Time-series (HOT) program: background, rationale and field implementation. Deep Sea Res. 2 43, 129–156 (1996).

    Google Scholar 

  30. Martiny, A. C., Vrugt, J. A. & Lomas, M. W. Concentrations and ratios of particulate organic carbon, nitrogen, and phosphorus in the global ocean. Sci. Data 1, 140048 (2014).

    Google Scholar 

  31. Martiny, A. C., Vrugt, J. A. & Lomas, M. W. Data from: Concentrations and ratios of particulate organic carbon, nitrogen, and phosphorus in the global ocean. Dryad (2015).

  32. Schwab, V. F. et al. 14C-free carbon Is a major contributor to cellular biomass in geochemically distinct groundwater of shallow sedimentary bedrock aquifers. Water Resour. Res. 55, 2104–2121 (2019).

    Google Scholar 

  33. Taubert, M. et al. Bolstering fitness via CO2 fixation and organic carbon uptake: mixotrophs in modern groundwater. ISME J 16, 1153–1162 (2022).

    Google Scholar 

  34. Rimstidt, J. D. & Vaughan, D. J. Pyrite oxidation: a state-of-the-art assessment of the reaction mechanism. Geochim. Cosmochim. Acta 67, 873–880 (2003).

    Google Scholar 

  35. Lin, W. et al. Genomic insights into the uncultured genus “Candidatus Magnetobacterium” in the phylum Nitrospirae. ISME J. 8, 2463–2477 (2014).

    Google Scholar 

  36. Kato, S. et al. Genome-enabled metabolic reconstruction of dominant chemosynthetic colonizers in deep-sea massive sulfide deposits. Environ. Microbiol. 20, 862–877 (2018).

    Google Scholar 

  37. Anantharaman, K. et al. Thousands of microbial genomes shed light on interconnected biogeochemical processes in an aquifer system. Nat. Commun. 7, 13219 (2016).

    Google Scholar 

  38. Kojima, H., Watanabe, T. & Fukui, M. Sulfuricaulis limicola gen. nov., sp. nov., a sulfur oxidizer isolated from a lake. Int. J. Syst. Evol. Microbiol. 66, 266–270 (2016).

    Google Scholar 

  39. Strous, M., Van Gerven, E., Kuenen, J. G. & Jetten, M. Effects of aerobic and microaerobic conditions on anaerobic ammonium-oxidizing (anammox) sludge. Appl. Environ. Microbiol. 63, 2446–2448 (1997).

    Google Scholar 

  40. Ji, X., Wu, Z., Sung, S. & Lee, P.-H. Metagenomics and metatranscriptomics analyses reveal oxygen detoxification and mixotrophic potentials of an enriched anammox culture in a continuous stirred-tank reactor. Water Res. 166, 115039 (2019).

    Google Scholar 

  41. Dalsgaard, T. et al. Oxygen at nanomolar levels reversibly suppresses process rates and gene expression in anammox and denitrification in the oxygen minimum zone off northern Chile. mBio 5, e01966 (2014).

    Google Scholar 

  42. Smith, R. L., Böhlke, J. K., Song, B. & Tobias, C. R. Role of anaerobic ammonium oxidation (anammox) in nitrogen removal from a freshwater aquifer. Environ. Sci. Technol. 49, 12169–12177 (2015).

    Google Scholar 

  43. Strous, M., Heijnen, J. J., Kuenen, J. G. & Jetten, M. S. M. The sequencing batch reactor as a powerful tool for the study of slowly growing anaerobic ammonium-oxidizing microorganisms. Appl. Microbiol. Biotechnol. 50, 589–596 (1998).

    Google Scholar 

  44. Kits, K. D. et al. Kinetic analysis of a complete nitrifier reveals an oligotrophic lifestyle. Nature 549, 269–272 (2017).

    Google Scholar 

  45. Rittmann, B. E. & McCarty, P. L. Environmental Biotechnology: Principles and Applications (McGraw-Hill Education, 2001).

  46. Zhang, Y. et al. Nitrifier adaptation to low energy flux controls inventory of reduced nitrogen in the dark ocean. Proc. Natl. Acad. Sci. USA 117, 4823–4830 (2020).

    Google Scholar 

  47. Field, C. B., Behrenfeld, M. J., Randerson, J. T. & Falkowski, P. Primary production of the biosphere: integrating terrestrial and oceanic components. Science 281, 237–240 (1998).

    Google Scholar 

  48. Lehmann, R. & Totsche, K. U. Multi-directional flow dynamics shape groundwater quality in sloping bedrock strata. J. Hydrol. 580, 124291 (2020).

    Google Scholar 

  49. Küsel, K. et al. How deep can surface signals be traced in the Critical Zone? Merging biodiversity with biogeochemistry research in a central German Muschelkalk landscape. Front. Earth Sci. 4, 32 (2016).

    Google Scholar 

  50. Yan, L. et al. Environmental selection shapes the formation of near-surface groundwater microbiomes. Water Res. 170, 115341 (2019).

    Google Scholar 

  51. Pack, M. A. et al. A method for measuring methane oxidation rates using low levels of 14C-labeled methane and accelerator mass spectrometry: methane oxidation rates by AMS. Limnol. Oceanogr. Methods 9, 245–260 (2011).

    Google Scholar 

  52. Nielsen, E. S. The use of radio-active carbon (C14) for measuring organic production in the sea. ICES J. Mar. Sci. 18, 117–140 (1952).

    Google Scholar 

  53. Xu, X. et al. Modifying a sealed tube zinc reduction method for preparation of AMS graphite targets: reducing background and attaining high precision. Nucl. Instrum. Methods Phys. Res. B 259, 320–329 (2007).

    Google Scholar 

  54. Merser, S. Acetabulum online interactive statistical calculators. Accessed Feb, 2021.

  55. Bermuda Oceanographic Timeseries, accessed 21 Oct 2020,

  56. Hawaiian Oceanographic Timeseries, accessed 21 Oct 2020,

  57. Hawaiian Oceanographic Timeseries, accessed 21 Oct 2020,

  58. Kumar, S. et al. Nitrogen loss from pristine carbonate-rock aquifers of the Hainich Critical Zone Exploratory (Germany) is primarily driven by chemolithoautotrophic anammox processes. Front. Microbiol. 8, 1951 (2017).

    Google Scholar 

  59. Füssel, J. et al. Nitrite oxidation in the Namibian oxygen minimum zone. ISME J. 6, 1200–1209 (2012).

    Google Scholar 

  60. McIlvin, M. R. & Altabet, M. A. Chemical conversion of nitrate and nitrite to nitrous oxide for nitrogen and oxygen isotopic analysis in freshwater and seawater. Anal. Chem. 77, 5589–5595 (2005).

    Google Scholar 

  61. Dalsgaard, T., Thamdrup, B., Farías, L. & Revsbech, N. P. Anammox and denitrification in the oxygen minimum zone of the eastern South Pacific. Limnol. Oceanogr. 57, 1331–1346 (2012).

    Google Scholar 

  62. Thamdrup, B. et al. Anaerobic ammonium oxidation in the oxygen-deficient waters off northern Chile. Limnol. Oceanogr. 51, 2145–2156 (2006).

    Google Scholar 

  63. Taubert, M. et al. Tracking active groundwater microbes with D2O labelling to understand their ecosystem function. Environ. Microbiol. 20, 369–384 (2018).

    Google Scholar 

  64. Bushnell, B. BBMap (SourceForge, 2014);

  65. Bornemann, T. L. V. et al. Geological degassing enhances microbial metabolism in the continental subsurface. Preprint at bioRxiv (2020).

  66. Nurk, S., Meleshko, D., Korobeynikov, A. & Pevzner, P. A. metaSPAdes: a new versatile metagenomic assembler. Genome Res. 27, 824–834 (2017).

    Google Scholar 

  67. Hyatt, D. et al. Prodigal: prokaryotic gene recognition and translation initiation site identification. BMC Bioinformatics 11, 119 (2010).

    Google Scholar 

  68. Langmead, B. & Salzberg, S. L. Fast gapped-read alignment with Bowtie 2. Nat. Methods 9, 357 (2012).

    Google Scholar 

  69. Wu, Y.-W., Simmons, B. A. & Singer, S. W. MaxBin 2.0: an automated binning algorithm to recover genomes from multiple metagenomic datasets. Bioinformatics 32, 605–607 (2016).

    Google Scholar 

  70. Brown, C. T. et al. Unusual biology across a group comprising more than 15% of domain bacteria. Nature 523, 208–211 (2015).

    Google Scholar 

  71. Sieber, C. M. K. et al. Recovery of genomes from metagenomes via a dereplication, aggregation and scoring strategy. Nat. Microbiol. 3, 836–843 (2018).

    Google Scholar 

  72. Olm, M. R., Brown, C. T., Brooks, B. & Banfield, J. F. dRep: a tool for fast and accurate genomic comparisons that enables improved genome recovery from metagenomes through de-replication. ISME J. 11, 2864–2868 (2017).

    Google Scholar 

  73. Murat Eren, A. et al. Anvi’o: an advanced analysis and visualization platform for ‘omics data. PeerJ 3, e1319 (2015).

    Google Scholar 

  74. Aramaki, T. et al. KofamKOALA: KEGG Ortholog assignment based on profile HMM and adaptive score threshold. Bioinformatics 36, 2251–2252 (2020).

    Google Scholar 

  75. Graham, E. D., Heidelberg, J. F. & Tully, B. J. Potential for primary productivity in a globally-distributed bacterial phototroph. ISME J. 12, 1861–1866 (2018).

    Google Scholar 

  76. Kanehisa, M., Sato, Y. & Morishima, K. BlastKOALA and GhostKOALA: KEGG tools for functional characterization of genome and metagenome sequences. J. Mol. Biol. 428, 726–731 (2016).

    Google Scholar 

  77. Altschul, S. F., Gish, W., Miller, W., Myers, E. W. & Lipman, D. J. Basic local alignment search tool. J. Mol. Biol. 215, 403–410 (1990).

    Google Scholar 

  78. Pelikan, C. et al. Diversity analysis of sulfite- and sulfate-reducing microorganisms by multiplex dsrA and dsrB amplicon sequencing using new primers and mock community-optimized bioinformatics. Environ. Microbiol. 18, 2994–3009 (2016).

    Google Scholar 

  79. Lücker, S., Nowka, B., Rattei, T., Spieck, E. & Daims, H. The genome of Nitrospina gracilis Illuminates the metabolism and evolution of the major marine nitrite oxidizer. Front. Microbiol. 4, 27 (2013).

    Google Scholar 

  80. Orellana, L. H., Rodriguez-R, L. M. & Konstantinidis, K. T. ROCker: accurate detection and quantification of target genes in short-read metagenomic data sets by modeling sliding-window bitscores. Nucleic Acids Res. 45, e14 (2017).

    Google Scholar 

  81. Chaumeil, P.-A., Mussig, A. J., Hugenholtz, P. & Parks, D. H. GTDB-Tk: a toolkit to classify genomes with the Genome Taxonomy Database. Bioinformatics 36, 1925–1927 (2020).

    Google Scholar 

  82. Parks, D. H. et al. A standardized bacterial taxonomy based on genome phylogeny substantially revises the tree of life. Nat. Biotechnol. 36, 996–1004 (2018).

    Google Scholar 

  83. Parks, D. H. et al. A complete domain-to-species taxonomy for Bacteria and Archaea. Nat. Biotechnol. (2020).

  84. Matsen, F. A., Kodner, R. B. & Armbrust, E. V. pplacer: linear time maximum-likelihood and Bayesian phylogenetic placement of sequences onto a fixed reference tree. BMC Bioinformatics 11, 538 (2010).

    Google Scholar 

  85. Jain, C., Rodriguez-R, L. M., Phillippy, A. M., Konstantinidis, K. T. & Aluru, S. High throughput ANI analysis of 90 K prokaryotic genomes reveals clear species boundaries. Nat. Commun. 9, 5114 (2018).

    Google Scholar 

  86. Price, M. N., Dehal, P. S. & Arkin, A. P. FastTree 2—approximately maximum-likelihood trees for large alignments. PLoS ONE 5, e9490 (2010).

    Google Scholar 

  87. Eddy, S. R. Accelerated profile HMM searches. PLoS Comput. Biol. 7, e1002195 (2011).

    Google Scholar 

  88. Ondov, B. D. et al. Mash: fast genome and metagenome distance estimation using MinHash. Genome Biol. 17, 132 (2016).

    Google Scholar 

  89. Letunic, I. & Bork, P. Interactive Tree Of Life (iTOL): an online tool for phylogenetic tree display and annotation. Bioinformatics 23, 127–128 (2007).

    Google Scholar 

  90. Emiola, A. & Oh, J. High throughput in situ metagenomic measurement of bacterial replication at ultra-low sequencing coverage. Nat. Commun. 9, 4956 (2018).

    Google Scholar 

  91. Wegner, C.-E. et al. Biogeochemical regimes in shallow aquifers reflect the metabolic coupling of the elements nitrogen, sulfur, and carbon. Appl. Environ. Microbiol. 85, e02346-18 (2019).

    Google Scholar 

  92. R: A Language and Environment for Statistical Computing (R Core Team, 2018).

  93. RStudio: Integrated Development Environment for R (RStudio Team, 2016).

  94. Wickham, H. et al. Welcome to the tidyverse. J. Open Source Softw. 4, 1686 (2019).

    Google Scholar 

  95. Neuwirth, E. RColorBrewer: ColorBrewer Palettes. R package version 1.1-2. (2014).

Download references


We thank F. Gutmann, H. Minkmar, P. A. Figueroa‐Gonzalez and R. Lehmann for their assistance with sample preparation, collection and filtration. We also thank M. Nowak and J. Southon for advice regarding 14C-bicarbonate supplementation and quantification of AMS results, respectively. In addition, we thank I. Görlich and M. Groth from the Core Facility DNA sequencing of the Leibniz Institute on Aging—Fritz Lipmann Institute in Jena for their help with Illumina sequencing. This study is part of the Collaborative Research Centre AquaDiva of the Friedrich Schiller University Jena, funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation)—SFB 1076—project number 218627073. S.T., X.X. and V.F.S. acknowledge additional support from the European Research Council (Horizon 2020 Research and Innovation Programme, grant agreement 695101). M.T. gratefully acknowledges funding from the DFG under Germany’s Excellence Strategy—EXC 2051—project ID 390713860. T.L.V.B. and A.J.P. were supported by the Ministerium für Kultur und Wissenschaft des Landes Nordrhein-Westfalen (‘Nachwuchsgruppe Dr. Alexander Probst’). Climate chambers to conduct experiments under controlled temperature conditions were financially supported by the Thüringer Ministerium für Wirtschaft, Wissenschaft und Digitale Gesellschaft (TMWWDG; project B 715-09075). The scientific results have in part been computed at the High-Performance Computing (HPC) Cluster EVE, a joint effort of both the Helmholtz Centre for Environmental Research—UFZ ( and the German Centre for Integrative Biodiversity Research (iDiv) Halle-Jena-Leipzig (


Open access funding provided by Friedrich-Schiller-Universität Jena.

Author information

Authors and Affiliations



K.K., S.T., W.A.O. and M.T. designed this study. W.A.O., K.K., A.J.P., T.L.V.B., M.M., M. Hölzer and M.T. planned, designed and conducted the metagenomic sampling approach. W.A.O., T.L.V.B. and A.J.P. performed the metagenomic analysis. W.A.O., S.T., X.X. and V.F.S. performed the carbon fixation experiments and interpreted the results. M.K., M. Herrmann, B.T. and L.B. conducted the anammox and nitrification rate measurements. W.A.O., K.K. and S.T. wrote the manuscript with the help of all authors.

Corresponding author

Correspondence to Kirsten Küsel.

Ethics declarations

Competing interests

The authors declare no competing interests.

Peer review

Peer review information

Nature Geoscience thanks Michael Wilkins and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Primary Handling Editors: Kyle Frischkorn and Tamara Goldin, in collaboration with the Nature Geoscience team.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended Data Fig. 1 Relative abundances of putative chemolithoautotrophic MAGs.

(Top) Proportion of metagenomic short reads that mapped to the putatively chemolithoautotrophic MAGs. Colors indicate the phylum and each bar indicates a specific MAG. (Bottom) The normalized coverages assigned to the MAGs, which unlike the proportional representation is independent of MAG size. For each plot, 0.1 µm and 0.2 µm represent the filter fraction from which the DNA was extracted from.

Extended Data Fig. 2 Relative abundances of transcripts mapping to putative chemolithoautotrophic MAGs.

Top panels are reproduced from Fig. 1 B, C. Middle panels are the corresponding values for the RNA-based reads. These values are derived from the ORF-based average normalized transcriptional activity of all putatively chemolithoautotrophic MAGs, as described in the methods and represent the sample-wise sum of the pie charts shown in Fig. 3. Bottom panel shows the average proportion of metatranscriptomic reads that map to the putatively chemolithoautotrophic MAGs. We did not have metatranscriptomes from H14 and the single H32 sample did not recruit reads.

Extended Data Fig. 3 Nitrification and anammox linear regressions.

Assessment of nitrogen transformation process rates in groundwater incubations by linear regression of measured increase in 15N over time for (a) nitrification with 15NH4+ label, (b) anammox with 15NH4+ label and (c) anammox with 15NO2- label. Shown are the slopes of the regression line for each replicate and control sample from which nitrification and anammox rates were calculated. R2 values describe the accuracy of the linearity and p-values the significance of results for p < 0.05.

Extended Data Fig. 4 Transcriptional activity of Sulfurifustaceae MAGs.

Transcriptional activity of both dominant and putatively chemolithoautotrophic MAGs within the family Sulfurifustaceae. These bins exhibited nearly identical distributions throughout the aquifer transect and a > 94% AAI. H32-bin014 is also shown in Fig. 4b. Bubble size is relative to the normalized ORF transcriptional coverage.

Extended Data Fig. 5 Phylogenetic tree of putative anammox MAGs.

Approximately maximum-likelihood phylogenetic tree of Wood Ljungdahl containing anammox MAGs and reference genomes. The multiple protein sequence alignment from GTDB_TK (5040 positions) was used with the JTT + CAT model, manually rooted with a MAG classified to the sister class Phycisphaerae within the phylum Planctomycetota.

Extended Data Table 1 CO2 fixation rates
Extended Data Table 2 Predicted global rates
Extended Data Table 3 Normalized rates

Supplementary information

Supplementary Information

Supplementary results and discussion.

Supplementary Table 1

MAG statistics, relative abundance profiles and transcriptional activity estimates

Supplementary Table 2

Raw and summarized AMS data

Supplementary Table 3

Sample information for all metagenomes constructed

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Overholt, W.A., Trumbore, S., Xu, X. et al. Carbon fixation rates in groundwater similar to those in oligotrophic marine systems. Nat. Geosci. 15, 561–567 (2022).

Download citation

  • Received:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI:


Quick links

Nature Briefing

Sign up for the Nature Briefing newsletter — what matters in science, free to your inbox daily.

Get the most important science stories of the day, free in your inbox. Sign up for Nature Briefing