Cold seeps are potential hotspots of deep-sea nitrogen loss driven by microorganisms across 21 phyla

Domain: earth_science
Published: 2025-02-14
Editorial outcome: minor_revision

Abstract
Nitrogen bioavailability, governed by fixation and loss processes, is crucial for oceanic productivity and global biogeochemical cycles. The key nitrogen loss organisms—denitrifiers and anaerobic ammonium-oxidizing (anammox) bacteria—remain poorly understood in deep-sea cold seeps. This study combined geochemical measurements, 15N stable isotope tracer analysis, metagenomics, metatranscriptomics, and three-dimensional protein structural simulations to explore cold-seeps nitrogen loss processes. Geochemical evidence from 359 sediment samples shows significantly higher nitrogen loss rates in cold seeps compared to typical deep-sea sediments, with nitrogen loss flux from surface sediments estimated at 4.96–7.63 Tg N yr-1 (1.65–2.54% of global marine sediment). Examination of 147 million non-redundant genes indicates a high prevalence of nitrogen loss genes, including nitrous-oxide reductase (NosZ; 6.88 genes per million reads, GPM), nitric oxide dismutase (Nod; 1.29 GPM), and hydrazine synthase (HzsA; 3.35 GPM) in surface sediments. Analysis of 3,164 metagenome-assembled genomes expands the nitrous-oxide reducers by three phyla, nitric oxide-dismutating organisms by one phylum and two orders, and anammox bacteria by ten phyla going beyond Planctomycetota. These microbes exhibit structural adaptations and complex gene cluster enabling survival in cold seeps. Cold seeps likely are previously underestimated nitrogen loss hotspots, potentially contributing notably to the global nitrogen cycle.

Introduction
Cold seeps are specialized marine environments primarily located along continental slopes and subduction zones, where subsurface fluids rich in hydrogen sulfide and hydrocarbons, such as methane, seep through the seabed1–3. These environments support complex ecosystems centered on the anaerobic oxidation of methane (AOM), a process involving methane-consuming archaea and sulfate-reducing bacteria. All life in these deep-sea oases depends on bioavailable nitrogen to support growth, which is a critical factor limiting biological productivity. Thus, understanding the processes that balance the nitrogen budget in cold seeps is essential4,5. Diazotrophs, organisms that convert dinitrogen gas (N2) into bioavailable ammonium (NH4+) through biological nitrogen fixation, are widespread in these environments, supported by diverse energy sources from either cultivated or uncultivated lineages6. Concurrently, nitrogen loss microbes convert bioavailable nitrogen back into N2 to maintain a balanced nitrogen cycle. However, studies on nitrogen loss processes in deep-sea cold seeps and the responsible microbial communities are relatively limited7–9. Understanding these processes is crucial for comprehending how cold seeps contribute to the broader nitrogen cycle in marine systems.

In marine sediments, two primary microbial processes for nitrogen loss are denitrification and anaerobic ammonium oxidation (anammox)10,11. The relative contributions of denitrification and anammox to nitrogen loss vary across different marine sediments. Most studies have focused on areas within 1000 meters of water depth, where denitrification typically accounts for over 80% of nitrogen loss10,12,13. Meanwhile, in deeper sea sediments (>1000 meters water depths), anammox emerges as a crucial process, which contributes up to 50% to the N2 production rates measured in sediments of the Cascadia Basin14, Peru margin15, and East China Sea16. These observations indicate that both denitrification and anammox may play a role in nitrogen loss at cold seeps. However, the exact contributions of these processes in balancing the nitrogen budget at cold seep sites remain unclear.

Denitrification occurs in two forms: classical and oxygenic. Classical denitrification reduces NO3− to NO2− and then nitric oxide (NO), which is subsequently converted to nitrous oxide (N2O). This N2O is further reduced to N2 by the enzyme nitrous-oxide reductase (NOR)17. NOR, encoded by nosZ clade I and II, is the only enzyme that biologically converts N2O to N218. NosZ clade I is primarily found in some members of Alpha-, Beta-, and Gammaproteobacteria that possess the complete denitrification pathway. In contrast, NosZ clade II has been found in diverse bacterial groups that lack nitrite reductase genes (nirS or nirK), including Gemmatimonadetes, Verrucomicrobia, Gammaproteobacteria, Campylobacterota, Myxococcota, Planctomycetota and Chloroflexi19–22. While NosZ clade I is well-documented, recent studies suggest that NosZ clade II may play a more important role in reducing N2O in certain ecosystems than previously assumed19,21,23,24. Oxygenic denitrification simplifies this process by converting NO into N2 and molecular oxygen (O2) through the action of nitric oxide dismutase (Nod), bypassing the production of N2O25,26. The nod genes were first identified in the genomes of Candidatus Methylomirabilis oxyfera-like bacteria (also known as NC10 or nitrite-dependent methane-oxidizing bacteria)26–28, which use Nod to generate N2 and O2. The generated O2 is utilized to catalyze methane oxidation, there by coupling oxygenic denitrification with aerobic methanotrophic pathways in anoxic environments27,29,30. These organisms have been detected in globally distributed cold seeps and other seafloor sediments through metagenomic and amplicon sequencing of pmoA or 16S rRNA genes31–33. Other bacteria like the Gammaproteobacterium strain HdN1, and species from the genera Sediminibacterium and Algoriphagus (within the phylum Bacteroidota)34, which also possess nod genes, may produce oxygen via dismutation, suggesting the potential involvement of microbes in oxygen production at cold seeps.

Compared to denitrification, anammox combines NO2− with NH4+ to form N2 using carbon dioxide (CO2) or carbonate as the carbon source, providing another nitrogen removal pathway with lower organic carbon requirements10,35. This process includes the reduction of NO3− to NO2−, which then reacts with NH4+ to create hydrazine (N2H4), a highly reactive and toxic compound with a low redox potential which then is oxidized to N2. The hydrazine-forming reaction, facilitated by hydrazine synthase (Hzs) is essential for the anammox pathway due to its unique biochemical properties. Known anammox bacteria are all affiliated with the Planctomycetota phylum, specifically within five families in the Brocadiales order: Candidatus Scalinduaceae, Candidatus Brocadiaceae, Candidatus Anammoxibacteraceae, Candidatus Bathyanammoxibiaceae, and Candidatus Subterrananammoxibiaceae36–38. Although these lineages are widely distributed across marine ecosystems and found in sediments from diverse marine environments, pure cultures have not yet been obtained39–41. Further studies using amplicon sequencing of 16S rRNA, hydrazine dehydrogenase (hzo), and hzsB genes have discovered diverse anammox bacteria in deep-sea cold seep sediments of the Okhotsk Sea and the South China Sea8,9. These findings suggest that the deep-sea cold seep environment might harbor previously undescribed anammox bacteria outside the Planctomycetota phylum.

In this study, we aim to investigate the contributions of denitrifying and anammox microbial communities in nitrogen loss at cold seep habitats, along with their diversity. We first provide geochemical evidence for nitrogen losses in cold seeps, based on data from 324 sediment samples collected from three sites—Lingshui, Haima, and Site F (Fig. 1 and Supplementary Fig. 1). This evidence is further supported by measurements of denitrification and anammox activities, conducted through slurry incubation experiments with 15N-labeled tracers. Subsequently, we explore the genes associated with nitrogen losses (nosZ, nod, and hzsA) and the diversity of microbes linked to these processes. This is achieved through sequence- and structure-based analyses using a detailed gene and genome catalog compiled from 165 metagenomes from 16 cold seep sites. Our findings reveal that cold seeps are overlooked areas for nitrogen loss in marine sediments under high pressure and low temperatures. Nitrogen loss in this habitat is mediated by diverse microbial populations, including newly identified phyla of anammox bacteria, contributing to this process.

Results
Geochemical evidence for nitrogen loss in cold seeps
Our geochemical measurements provide evidence for nitrogen loss in cold seep sediments. A total of 324 sediment samples from 37 cores with depths of 0–36 cm below the seafloor (cmbsf) were collected from three different cold seeps—Lingshui, Haima, and Site F (Supplementary Fig. 1). These sites represent different stages of cold seep activity with varying methane seeping intensity42. Correspondingly, sulfate (SO42−) concentrations in the porewater showed coherent downcore variations, either decreasing rapidly or remaining relatively stable depending on methane fluxes (Fig. 1a and Supplementary Data 1). The predominant form of dissolved inorganic nitrogen was NH4+ (Fig. 1b), with the highest concentrations found at Site F (8.6–867.2 μM, average 266.0 μM) and lower concentrations at Lingshui and Haima (averaged at 61.5 and 15.0 μM, respectively). Porewater NO3− and NO2− concentrations (averaged at 1.6 and 0.6 μM, respectively) were generally lower than NH4+ across all sites (Fig. 1c-d). Following this, we detected the signature of active N2 production, as evident by the downcore increasing concentrations of N2 (1.89–6.12% of the headspace gas) in two cores, with δ15N values in the range of −3.08–6.11‰ (Fig. 1e, f). The depleted values of δ15N (−3.09–0.10‰) indicated intense denitrification and anammox activities at Site F-1443,44. Together, our extensive geochemical measurements indicate the occurrence of nitrogen loss in cold seep sediments.

Total organic carbon (TOC), total nitrogen (TN), and TOC:TN ratio varied considerably between sites, with both TOC and TN generally decreasing with sediment depth (Fig. 1g–i). TOC levels at the Lingshui-6/-10/-11/-12, Haima-7/-8/-9, and Site F-14/-15 sites were much higher than other sites, suggesting that these sites could provide more electron donors for heterotrophic denitrification and leading to potentially higher nitrogen loss rates (Supplementary Fig. 2). Enzymes involved in microbial nitrogen loss processes require metal cofactors45–47, which were abundant in all studied sediments. In Lingshui, copper (Cu2+) reached the highest concentration, averaging 26.6 nM, while zinc (Zn2+) peaked in Haima, averaging 100.4 nM (Fig. 1j, k). Calcium (Ca2+) displayed a decreasing trend in Site F, while iron (Fe2+) was highest in Haima, averaging 11.9 μM (Fig. 1l, m).

Slurry experiments indicate substantial nitrogen loss rates in cold seep sediments
We examined potential nitrogen loss rates and the contribution of anammox to N2 production by conducting slurry experiments with nitrogen isotope tracing on sediment samples up to 36 cm deep from Lingshui, Haima, Site F, Shenhu, and Qiongdongnan (Supplementary Fig. 1) at 4 °C. We detected both denitrification and anammox, with rates significantly higher in cold seep regions (n=60) compared to non-seep regions (n=34; P<0.0001; Fig. 2b and Supplementary Data 2). In cold seeps, average potential denitrification rates were measured at 2.30±4.03 nmol cm−3 h−1, and anammox rates were 0.14±0.21 nmol cm−3 h−1. In contrast, non-seep regions had denitrification and anammox rates of only 0.1±0.1 nmol cm−3 h−1 and 0.04±0.05 nmol cm−3 h−1, respectively. It should be noted that the lesser importance of anammox in the cold seep sediments is mainly attributed to higher organic carbon contents (Fig. 1g), which favor denitrification over anammox11,48.

Although the anammox contributions were lower than those measured previously in deep-sea sediments (Fig. 2a), they were comparable to values measured in sediments from continental slope15,41 and deep-sea14 environments. Denitrification is the primary nitrogen loss process in these environments, with anammox contributing to roughly 26.14±28.93% of the total nitrogen loss in cold seeps (Supplementary Fig. 2). The rich organic carbon and tight correlation between denitrification rates and TOC contents in cold seeps (P<0.01; Supplementary Fig. 3a) suggest that heterotrophic denitrification is a likely metabolic pathway7. However, deep-sea cold seeps are chemosynthetic ecosystems where seepage fluids rich in hydrogen sulfide, methane, and other hydrocarbons provide energy for diverse microbial communities49. These chemosynthetic communities may also play an important role in the denitrification process. For instance, chemosynthetic Campylobacterota are abundant in cold seep sediments (Fig. 3a), which are likely responsible for sulfide oxidation coupled with NO3− reduction50. In addition, Beggiatoaceae, a chemolithoautotrophic sulfur-oxidizing bacterium51, is one of the most prevalent NosZ-encoding organisms in cold seep sediments (Supplementary Data 3). A tight correlation between anammox and denitrification rates was observed in cold seeps and non-seep (P<0.01; Supplementary Fig. 3c and d), aligning with the fact that denitrification from nitrate to nitrite provides the necessary NO2− for anammox in marine sediments52,53.

Despite low-temperature conditions, several cold seep sites—specifically Haima-6, Haima-7, Haima-8, Lingshui-10, and Lingshui-11—exhibited exceptionally high rates of nitrogen loss (up to 17.65±1.24 nmol cm−3 h−1; Supplementary Fig. 2 and Supplementary Data 2) possibly related to abundant carbon sources (1.01–2.37%; Supplementary Data 1). These rates are comparable to those in estuarine and coastal environments of higher temperatures, as well as other known cold seep sites such as those in the Gulf of Mexico7, but are considerably higher than typical deep-sea sediments (Fig. 2a and Supplementary Data 3). However, it should be noted that these comparisons are approximations, as variations in organic matter concentrations and composition across regions introduce uncertainties that limit the reliability of direct comparisons.

These findings suggest that certain denitrifying or anammox microbes have adapted or evolved to thrive in deep-sea cold seeps with high methane and sulfide concentrations. Although our experiments did not replicate the high hydrostatic pressure found in natural settings, research has shown that N2 production rates measured in the laboratory closely matched those observed in situ using a benthic lander in hadal sediments of the Atacama Trench11. We extrapolated the potential nitrogen loss rates measured in our study across the globally estimated seepage area using methods from a previous study54. Assuming a sediment mean bulk dry density of about 1.3 g cm−355 and an estimated global cold seep area of 2.05–3.15 × 105 km256,57, we calculated a global nitrogen loss flux from denitrification and anammox in the top 5 cm of cold seep sediments to be around 4.96–7.63 Tg N per year. This flux represents roughly 1.65–2.54% of the global marine sediment nitrogen loss (estimated at 300 Tg N per year)58, underscoring that, despite covering only about 0.057–0.087% of the global marine area, cold seeps are potential nitrogen loss hotspots in deep-sea sediments. It is important to note that the measured rates are potential rates, meaning that they could be higher than actual in situ rates by a variable and unknown degree. Our use of slurry incubations provided a controlled method for assessing nitrogen loss potential but did not fully replicate natural conditions and, therefore, might not capture exact in situ rates. In addition, the extent of global cold seep areas remains imprecise, as cold seep sizes can vary significantly, adding uncertainty to our estimate. Nonetheless, this calculation provides a first-order quantification that highlights the potential importance of cold seeps in global sedimentary nitrogen loss. Future studies with in situ techniques and advanced modeling approaches will help refine these estimates and better capture the true nitrogen loss rates in these unique environments.

Near-surface sediments are hotspots of nitrogen loss gene diversity
To investigate the diversity of nitrogen loss genes in cold seeps, we further focus on their evolution, taxonomy, and potential horizontal transfer via mobile genetic elements (MGEs). Using a gene catalog of 147 million non-redundant genes from cold seeps59, we delved into the diversity of genes linked to nitrogen loss, focusing on nitrous-oxide reductase (NosZ), nitric oxide dismutase (Nod), and hydrazine synthase (HzsA) (Supplementary Fig. 4). We identified 530 NosZ sequences containing cupredoxin-related protein active domains (Supplementary Fig. 5). These sequences are categorized into two groups based on their signal peptides: clade I (n=164) utilizes the twin-arginine translocation (Tat) pathway60, while clade II (n=366) has an additional c-type heme domain at the C-terminus, associated with the secretory (Sec) pathway20,24. Furthermore, we discovered a sequence-divergent branch called NosZG7 (n=27), which, despite sequence differences, showed structural congruence with canonical NosZ (Supplementary Fig. 5b). In addition, we identified 151 putative Nod sequences, all bearing conserved catalytic site residues similar to those of the putative nitric oxide dismutase34 from the NO-dismutating bacterium Methylomirabilis oxyfera (Supplementary Fig. 6). The sequence-divergent branches Cluster1 (n=28) and Cluster5 (n=8) were only found through phylogenetic and structural analysis, respectively. We also identified 644 HzsA sequences with hydrazine synthase alpha subunit domains and the pentacoordinated c-type heme47. The HzsA sequences were classified into six clades, with five branches identified through structural analysis as divergent from those of known anammox bacteria (Supplementary Fig. 7). Among the genes of the nosZ clade II, the most abundant contig belongs to Gammaproteobacteria, while that of nosZ clade I genes is affiliated with Alphaproteobacteria (Supplementary Data 4). Furthermore, the predominant contigs of nod and hzsA genes are assigned to Scalindua (Planctomycetota) and Bacteroidales (Bacteroidota), respectively (Supplementary Data 4). Notably, while no hzsA and nod protein sequences were detected in MGEs, we found 31 nosZ sequences in MGEs, with clade I (n=24; Supplementary Fig. 8) being more prevalent than clade II (n=7). This observation suggests that horizontal gene transfer via MGEs may play a significant role in the diversification of nosZ-bearing denitrifying microorganisms in deep-sea cold seeps, supporting the hypothesis that MGEs are crucial for genetic exchange and adaptation in these environments61.

The abundance, expression, and depth distribution of nitrogen loss genes in cold seep sediments were explored, revealing their role in nitrogen cycling dynamics and their relationship with environmental drivers. The average abundances of nosZ, nod, and hzsA genes were 4.84, 0.90, and 2.78 genes per million reads (GPM; 0–6855 cmbsf; Fig. 2c, Supplementary Fig. 9a and Supplementary Data 5), respectively. These abundances are one-quarter that of the reductive dsrA gene (averaging 21.32 GPM), indicative of sulfate reduction, and only one-tenth of the oxidative mcrA gene (averaging 46.15 GPM), which is indicative of methane oxidation. Similarly, they are one-tenth the abundance of the nifH gene (averaging 55.45 GPM), associated with nitrogen fixation. However, our measurements found that sediment nitrogen fixation rates were significantly lower than nitrogen loss rates in cold seeps (P<0.001; Supplementary Fig. 10 and Supplementary Data 2). This apparent discrepancy can be explained by the energy-intensive nature of nitrogen fixation, which requires substantial energy and reducing equivalents62, thereby limiting its activity despite higher gene abundance. In addition, alternative nitrogen sources, like ammonium from subsurface seepage56, and efficient recycling processes such as ammonification43, sustain nitrogen availability, offsetting losses and maintaining productivity in cold seep ecosystems.

The expression levels of nitrogen loss genes indicate that microbial nitrogen loss processes are active, particularly at the surface of cold seep sediments, with average values of 2.38 transcripts per million reads (TPM) for nosZ, 0.20 TPM for nod, and 0.48 TPM for hzsA (Fig. 2d, Supplementary Fig. 9b and Supplementary Data 6). Specifically, nosZ clade II genes were more abundant than clade I genes (P<0.0001; Fig. 2c, Supplementary Fig. 9a and Supplementary Data 5). The expression level of nosZ clade II genes (averaging 3.71 TPM) was also higher than that of nosZ clade I genes (averaging 1.13 TPM; P<0.0001; Supplementary Fig. 9b), indicating that nosZ clade II may play a more important role in N2O consumption in cold seeps19,23,24. The distributions of nosZ, nod, and hzsA genes varied with sediment depth (0–300 cmbsf), showing generally negative trends with increasing depth, indicating diminished nitrogen loss activity deeper in the sediments (Fig. 2e). However, within the shallow surface layers (up to 40 cmbsf), the abundance of nosZ, nod, and hzsA genes positively correlated with depth (P<0.05; Supplementary Fig. 9c and Supplementary Data 8). In addition, cold seep nitrogen loss gene abundance exhibited statistically significant differences across surface and middle/deep depths, decreasing with sediment depth, with lower abundances in deeper sediments. Specifically, gene abundances ranged from 1.29 to 6.88 GPM at the surface (0–50 cmbsf), 0.08 to 1.41 GPM in the middle depth (50–500 cmbsf), and 0.20 to 0.98 GPM in deep sediments (>500 cmbsf) (P<0.0001; Fig. 2c, e and Supplementary Data 5). To identify environmental drivers in the shallow surface layers (up to 40 cmbsf), we correlated the relative abundance of nitrogen loss genes with environmental factors (Supplementary Fig. 9c). Overall, NH4+ and NO3− were the strongest correlates of nosZ, nod, and hzsA genes in the surface layer (P<0.05). The correlation points towards a dominant contribution of denitrifiers and anammox bacteria, which is linked to NH4+ and NO3− conversion.

NosZ and Nod proteins from seep denitrifiers exhibit considerable structural diversity
To explore the diversity, distribution, and structural features of nosZ genes in cold seep environments, we then focus on their phylogenetic variation, transcriptional activity, and metabolic diversity across different microbial phyla. From a cold seep genome catalog of 3164 metagenome-assembled genomes (MAGs), we identified 142 nosZ sequences, which were categorized into nosZ clade I (n=52) and nosZ clade II (n=90) based on phylogenetic and structural analyses (Supplementary Fig. 11 and Supplementary Data 9). The relative abundance of the MAGs that comprised nosZ genes ranged from 0.0001 to 0.0263% in the 165 sediment samples of cold seeps (Supplementary Data 9). Structurally, the N-terminal of cold seep nosZ clades I and II NORs configurations differ to align with their physiological functions (Fig. 3a), adopting signal peptides63 for either the Tat pathway in clade I or the Sec pathway in clade II (Supplementary Figs. 12, 13 and Supplementary Data 15). In addition, some nosZ clade II variants possess a C-terminal α-helix, which surrounds the periphery of CuA and CuZ sites of adjacent monomers in the structure (Fig. 3a and Supplementary Fig. 13b), potentially enhancing the stability of their active sites. Details of protein folding are provided in the Supplementary Note 1.

The distribution of the 142 nosZ genes across MAGs spans one archaeal and 18 bacterial phyla, reflecting the wide phylogenetic breadth of potential nitrous-oxide reducers in these environments (Supplementary Fig. 11a and Supplementary Data 9). The most common phyla containing nosZ genes are Pseudomonadota (n=56) and Bacteroidota (n=36), with fewer occurrences in Campylobacterota, Myxococcota, Chloroflexota, Desulfobacterota, Gemmatimonadota and other phyla (Fig. 3a). NosZ clade I genes are present in Alphaproteobacteria, Gammaproteobacteria, and a few other phyla, while clade II genes are more widespread, found in 15 bacterial groups and the archaeal phylum Thermoplasmatota (Supplementary Fig. 11 and Supplementary Data 9). Among the MAGs containing nosZ clade II genes, the most abundant MAG belongs to Desulfobacterota (Desulfobacterales, 0.0263%), whereas the most abundant MAG with nosZ clade I genes is affiliated with Pseudomonadota (Photobacterium frigidiphilum, 0.0137%). Desulfobacterales and Photobacterium frigidiphilum function as chemoautotrophs64 and heterotrophs65,66, respectively, highlighting the metabolic diversity of denitrifying microbes in cold seeps. This reveals three additional phyla (Desulfobacterota, Krumholzibacteriota, and Zixibacteria) capable of reducing N2O, considerably expanding the known genetic diversity of potential N2O reducers in deep-sea cold seep sediments. Most clade II MAGs contain multiple nosZ genes, with six MAGs containing two copies and one containing three copies, contrasting with the single nosZ gene typically found in clade I MAGs from cold seeps (Fig. 3b, Supplementary Fig. 14 and Supplementary Data 13). Notably, we identified a MAG, RS_10_sbin_88, from Gammaproteobacteria that carries one copy of both nosZ clade I and II genes. The nosZ clade I genes of Pseudomonadota were transcribed at high levels, up to 150.71 TPM, while nosZ clade II genes were transcribed at moderate to high levels in several phyla, up to 41.92–200.31 TPM, especially in Campylobacterota and Pseudomonadota (Supplementary Data 10, 11).

In addition, the diversity, distribution, and structure of nod genes involved in nitric oxide dismutation in cold seeps were explored. From 3164 MAGs, we identified five nod genes characterized by conserved enzymatic active center residues, distinguishing them from the related nitric oxide reductase (nor) by the substitution of “Thr” with “Ile”, “His” with “Asp”, and “Glu” with “Gln”34 (Fig. 4). The relative abundance of MAGs comprising nod genes ranged from 0.0004 to 0.0133% in the 165 samples of cold seeps (Supplementary Data 16). Structurally, the nitric oxide dismutase from cold seeps consistently exhibited conserved α-helix domains—four enzymes had 14 α-helices, and one had 13, reflecting notable structural conservation (Supplementary Fig. 15a–e and Supplementary Data 17). These nod genes were found in five MAGs across two bacterial phyla, indicating a broader presence of oxygenic denitrifiers than previously understood (Fig. 4a, b and Supplementary Data 16)25,27,34. Specifically, three of these MAGs are attributed to the Planctomycetota phylum (two from UBA1135 and one from Scalinduaceae), and two to the Bacteroidota phylum (as Maribacter_A sp023141835 and Cecembia rubra from the orders Flavobacteriales and Cytophagales, respectively). This expanded diversity suggests more widespread bacterial capabilities for nitric oxide dismutation at cold seeps than previously appreciated through metagenomic or environmental studies67. The gene cluster for nitric oxide dismutase in UBA1135 is conserved across two genomes, featuring a helix-turn-helix domain and electron transport-associated proteins, including ferredoxin oxidoreductase and a 4Fe-4S dicluster domain (Supplementary Data 18). In Cecembia rubra, the nod gene cluster includes not only the standard nos cluster (nosZ, nosD, nosF, nosY, and nosL) but also elements that regulate nitric oxide signaling and electron transport. However, in the other two MAGs from Scalinduaceae and Maribacter_A sp023141835, the genes downstream of nod are categorized as having domains of unknown function. In addition, most of these putative nod-containing MAGs possess genes for other nitrogen metabolic processes (Supplementary Data 16) such as genes for hydroxylamine dehydrogenase (hao), nitric oxide reduction (nor), nitrous oxide reduction (nosZ), nitrate reduction (nap/nir), nitrite reduction to nitric oxide (nir) and nitrite reduction to ammonia (nrf/nir).

Genetic potential for anammox is found in multiple phyla beyond Planctomycetota
To explore the diversity, distribution, and structural characteristics of hzsA genes involved in anammox processes in cold seep sediments, further analyses of the 3,164 MAGs identified 265 hzsA genes, the diagnostic gene of bacteria capable of anammox (Fig. 5a and Supplementary Data 19). The relative abundances of these hzsA-contain MAGs ranged from 0.0002% to 0.0201% in the 165 samples from cold seeps (Supplementary Data 19). These hzsA genes sequences, characterized by a motif binding to a pentacoordinated c-type heme in the hydrazine synthase alpha subunit47, were distributed across 94 bacterial MAGs spanning 10 bacterial phyla (Supplementary Fig. 16 and Supplementary Data 19). This distribution extends beyond the previously recognized anammox bacteria, which contain traditional HzsA found in the Brocadiales order of the Planctomycetota phylum36–38,68, to include previously undescribed HzsA in species outside of Brocadiales. These hzsA-containing MAGs are members of Planctomycetota (n=55), Bacteroidota (n=16), Acidobacteriota (n=10), Verrucomicrobiota (n=6), Sumerlaeota (n=1), JABMQX01 (n=1), Calditrichota (n=1), Desulfobacterota (n=1), Fibrobacterota (n=1), Gemmatimonadota (n=1) and Zixibacteria (n=1) (Supplementary Fig. 16 and Supplementary Data 19). The hzsA genes of Acidobacteriota, Bacteroidota, Calditrichota, Planctomycetota, and Sumerlaeota were transcribed at moderate to high levels, up to 56.75–414.23 TPM, whereas fewer transcripts from Verrucomicrobiota were detected (Supplementary Data 20).

Among the 94 MAGs, 27 contain 42 hzd genes, which encode hydrazine dehydrogenase for oxidizing N2H4 to N2 (Supplementary Data 22). Notably, only two belong to the Ca. Scalinduaceae family within the Planctomycetota phylum, and both do not contain the complete hzsBC genes (Fig. 6 and Supplementary Data 23). However, 37 MAGs do contain these genes (Supplementary Data 23, 24), forming a multienzyme complex HZS-αβγ69. Of these, five MAGs comprise multiple hzsABC clusters and are associated with Phycisphaerae (n=3), Verrucomicrobiae (n=1), and Bacteroidia (n=1). These findings considerably broaden the recognized phylogenetic diversity and environmental presence of potential anammox bacteria.

We carefully examined the five MAGs with multiple hzsABC clusters and two Scalinduaceae MAGs, confirming they were not contaminated during binning (Supplementary Data 25). We found that the MAGs not only shared a similar gene cluster arrangement (Fig. 6a and Supplementary Data 23) but also structural similarities with known anammox bacteria (Fig. 5a). In addition to the essential genetic machinery for anammox metabolism, including hydrazine synthase (hzs) for converting NO and NH4+ to N2H4, and hydrazine dehydrogenase (hzd) for oxidizing N2H4 to N2, the MAGs also possess additional nitrogen metabolic processes such as genes for nitrate reduction (nap/nir), nitrite reduction to nitric oxide (nir) and nitrite reduction to ammonia (nrf/nir)70 (Fig. 6b and Supplementary Data 26). Notably, we identified a MAG affiliating with Verrucomicrobiae, FR_S1_sbin_24, which contained at least two copies of the hzsA, hzsB, and hzsC genes in a single contig (Fig. 6a and Supplementary Data 23), interspersed with several unknown proteins and upstream genes related to c-type cytochromes and Ca-activated chloride channels, which were annotated as belonging to the Verrucomicrobiae against NCBI’s non-redundant database. This MAG also encodes both hydroxylamine oxidoreductase (Hao) for nitrite reduction71 and NosZ (Fig. 6b), suggesting a versatile metabolic capability.

Both traditional and previously undescribed HzsA proteins encoded in the seven MAGs are structured into three distinct domains (Fig. 5a, Supplementary Figs. 17a, b, and 18a, b and Supplementary Note 2)47. However, traditional HzsA proteins lack an N-terminal transmembrane α-helix (Supplementary Fig. 17c, d), as do the previously undescribed HzsA proteins described here (Supplementary Fig. 18c, d). Only the previously undescribed HzsA proteins include N-terminal signal peptides indicative of the Sec secretion pathway (Supplementary Fig. 18e, f). In addition, these proteins exhibit significant variation in their C-terminal structures, particularly the previously undescribed HzsA proteins, which include a C-terminal α-helix that may stabilize the binding at the heme αI active center (Fig. 5a and Supplementary Figs. 17a, b, 18a, b). These proteins assemble into one or two heterotrimer complexes as predicted by AlphaFold Multimer (Fig. 5c and Supplementary Fig. 19). In each complex, one of the β and γ subunits of the previously undescribed HzsA are fused into a single polypeptide chain, consistent with earlier studies72. The previously undescribed HzsA proteins retain essential domains for heme binding and electron transfer, indicating they are functionally capable of catalyzing anammox reactions. Structural adaptations, such as N-terminal signal peptides and C-terminal α-helices, may optimize these proteins for the unique conditions of cold seeps, showcasing the evolutionary adaptability of anammox bacteria in nitrogen cycling.

Implications of nitrogen loss in cold seeps
In the past, gene annotations primarily relied on sequence similarity, which often overlooked genes that are distantly related yet functionally similar. However, over long evolutionary time scales, multiple substitutions at the same site can cause uncertainty in sequence alignment. Structures evolve at a slower rate than the underlying sequence mutations and are more conserved, emphasizing the importance of protein structural analysis73. By employing protein structural similarities and phylogenetic analysis, we have discovered that the NosZ clade II is notably more diverse and abundant in deep-sea cold seeps than previously assumed. This diversity may have significant implications for the role of NosZ clade II in N2O consumption in cold seeps. Members of the Planctomycetota phylum, as well as the orders Flavobacteriales and Cytophagales, might be notable contributors to nitrogen loss through nitric oxide dismutation. These organisms might also act as oxygen producers which could be linked to the presence of aerobic methane and sulfide oxidizers within apparently anoxic layers of cold seep ecosystems74. Crucially, our findings also indicate that anammox bacteria are found in multiple phyla (e.g., Bacteroidota, Acidobacteriota, and Verrucomicrobiota) beyond Planctomycetota and that these overlooked lineages actively express anammox genes. Our findings highlight the importance of nitrogen loss through denitrification and anammox processes at cold seeps. Overall, this study provides evidence supporting the presence of numerous previously undescribed, cold-adapted microbial lineages involved in denitrification (including oxygenic denitrification) and anammox processes. It establishes cold seeps as considerable nitrogen loss hotspots in the deep sea and as important contributors to the global nitrogen cycle, broadening the recognized roles of cold seeps beyond serving as oases of diversity, productivity, and methane removal.

Methods
Sampling and geochemical measurements
A total of 324 samples from 37 push cores were collected by a remotely operated vehicle (ROV) at Lingshui, Haima, and Site F cold seeps in the South China Sea between 2020 and 2024 (Supplementary Fig. 1). Upon retrieval on deck, sediment cores were immediately placed in a helium-filled glove bag. Porewater was extracted at 2 cm depth intervals using Rhizon samplers with 0.2 μm pore size (Rhizosphere, Netherlands). Collected porewater samples for use for metal analysis were acidified to pH ~2 with HNO3 (Optima grade, Thermo Scientific, USA), and stored at 4 °C before analysis. For nitrogen gas measurement in sediments, 3 mL of sediment at 2 cm intervals was transferred using a 5 mL cut-off syringe into a 22 mL serum vial. The vial was immediately crimp-sealed and allowed to equilibrate to enable the dissolved gases to partition into the headspace at 4 °C. Samples of porewater and sediment for analyzing nutrients, sulfate, and calcium were stored at −20 °C until analysis.

The concentrations of N2 in the headspace gas, along with its δ15N values were measured using gas chromatography with thermal conductivity detection (Agilent, USA) and a continuous-flow isotope-ratio mass spectrometer (SerCon, UK). The N2 concentration is expressed as the percentage of N2 in the total headspace gas volume of the 22 mL serum vial. The concentrations of ammonium (NH4+), nitrite (NO2−), and nitrate (NO3−) in porewater were determined using a Quattro continuous flow analyzer (SEAL Analytical AA3, Germany). Sulfate (SO42−) and calcium (Ca2+) concentrations were determined using a Dionex Ion Chromatograph (Thermo Scientific Dionex, USA). Dissolved metals (Fe2+, Cu2+, and Zn2+) in porewater were measured using ICP-MS (Thermo Scientific, USA). Sediment total organic carbon (TOC) and total nitrogen (TN) were measured using a Vario Micro Cube element analyzer (Elementar, Germany) after the sediments were treated with 1 M HCl to remove carbonates.

Determination of nitrogen loss rates
Nitrogen loss rates were determined using 60 samples from 12 push cores (Lingshui-10/-11/-12, Haima-3/-4/-5/-6/-7/-8/-9, and Site F-15/-16) collected by an ROV during 2023 and 2024 from the cold seeps of Lingshui, Haima, and Site F, with water depths ranging from 1171 to 1758 meters. In addition, 34 samples from 8 push cores (Shenhu-1/-2/-3/-4/-5, Qiongdong-1, Lingshui-13, and Site F-17) were collected from the non-seep area in Shenhu, Lingshui, Site F, and Qiongdongnan, at depths of 489 to 1761 meters in the South China Sea (Supplementary Fig. 1). Potential nitrogen loss rates were measured using N isotope-tracing techniques as follows75,76.

Potential denitrification and anammox rates were measured using the N-isotope tracing technique based on previous studies with some modifications48,75,77. Briefly, slurries were prepared using collected sediments and artificial seawater matching in situ salinity, at a sediment/water volume ratio of 1:7. The mixture was purged with helium for 30 min and stirred vigorously to ensure homogeneity. Gas-tight borosilicate vials (n=21; Exetainers, Labco, UK) were then filled with slurries under a helium atmosphere. Subsequently, the vials were pre-incubated at near in situ temperature (4 °C) for 48 h to remove residual NO3− (NO3−+NO2−) and dissolved oxygen. After the pre-incubation, three vials were used to measure the ambient NO3−. The residual slurries (n=18) were divided into three groups, each amended with different nitrogen compounds in helium-purged stock solutions: (1) 15NH4+ (99.12%), (2) 15NH4+ + 14NO3−, and (3) 15NO3− (99.21%). The final concentration of 15N compounds in each vial was 75 μM. At 0 and 8 h incubation times, triplicate vials were sacrificed with 200 μL 50% ZnCl2 solution as initial and final samples, respectively. The 15N labeled N2 were measured using a membrane inlet mass spectrometer (MIMS, HPR-20, Hiden Analytical, UK) with a detection limit of 0.1 μM. Denitrification and anammox rates were estimated from the accumulation of 29N2 and 30N2 during the slurry incubation75,78. The respective contributions of denitrification and anammox to 29N2 production were quantified using Eq. (1).

P29 = A29 + D29 (1)

where P29 (nmol cm−3 h−1) denotes the total 29N2 production rates, D29 (nmol cm−3 h−1) and A29 (nmol cm−3 h−1) denote the production rates of 29N2 from denitrification and anammox, respectively. Here, D29 was obtained by Eq. (2), assuming random pairing of 14N and 15N from 14NO3− or 15NO3−79,80.

D29 = P30 × 2 × (1 − FN) × FN − 1 (2)

where P30 (nmol cm−3 h−1) denotes the total 30N2 production rates, FN (%) denotes the fraction of 15N in NO3−, which was obtained from the added 15NO3− and the measured residual ambient NO3−, ranging from 97.8% to 99.7%. The potential rates of denitrification and anammox were quantified by Eqs. (3) and (4).

Dt = D29 + 2 × P30 (3)

At = P29 − D29 (4)

where Dt and At (nmol cm−3 h−1) denote the denitrification and anammox rates, respectively. By convention, the percent of N2 production accounted for anammox is abbreviated as ra (%).

Metagenomic data processing
The metagenomic datasets from 165 samples were collected from 16 cold seep sites worldwide, including oil and gas seeps, methane seeps, gas hydrates, asphalt volcanoes, and mud volcanoes (Supplementary Fig. 1). Non-redundant gene and genome catalogs were constructed in our previous study59. Briefly, metagenomic sequence data were quality-controlled and assembled into contigs. Protein-coding sequences were then predicted and clustered to create a non-redundant gene catalog consisting of 147,289,169 representative clusters. Salmon (v1.10.2)81 was used to calculate gene abundance in each metagenome, which was then normalized to genes per million reads (GPM)82. Functional annotations were performed using eggNOG-mapper (v2.1.9) with default parameters83,84. MMseqs2 taxonomy (v13.45111; parameter: --tax-lineage 1)85 was used to assign taxonomic labels to each representative amino acid sequence with reference to GTDB R207 database86.

Contigs longer than 1000 bp were selected for subsequent binning, and the produced metagenome-assembled genomes (MAGs) underwent dereplication at 95% average nucleotide identity59. Then, the MAGs were checked by GUNC (v1.0.5; default parameters)87 to remove genomes potentially containing chimerism, retaining only those classified as “pass.GUNC”. A total of 3164 representative MAGs (completeness >50% and contamination <10%) were obtained. The relative abundance of each MAG was calculated using CoverM (v0.6.1; https://github.com/wwood/CoverM; parameters: -m relative_abundance --trim-min 0.10 --trim-max 0.90 --min-read-percent-identity 0.95 --min-read-aligned-percent 0.75). The taxonomy of each MAG was assigned using GTDB-Tk v2.1.1 with reference to GTDB R207 database86.

Identification of nitrous oxide reductase gene (nosZ)
For nosZ gene database search (see workflow in Supplementary Fig. 4a), we first queried protein sequences from the non-redundant gene catalog and 3164 MAGs against NCycDB88 using DIAMOND (version 2.0.14)89 in blastp mode (-k 1 -e 0.0001 -p 5). Subsequently, the nosZ reference sequences (n=403) from the Greening lab metabolic marker gene databases90 were used to search for potential nosZ sequences in the non-redundant gene catalog and 3164 MAGs with DIAMOND blastp (version 2.0.14; --id 50)89. In addition, hidden Markov models (HMMs) of nosZ clade I (TAT-dependent nitrous-oxide reductase) and nosZ clade II (Sec-dependent nitrous-oxide reductase) were obtained from NCBI’s Protein Family Models using accession “TIGR04244.1”91 and “TIGR04246.1”92, respectively. These HMMs were used to screen proteins from the non-redundant gene catalog and 3164 MAGs with hmmsearch in HMMER v3.3.2 using the parameter -E 1e-5. All nosZ genes identified with the above three methods were merged, and any sequences shorter than 400 amino acids were excluded.

The phylogenetic trees of the amino acid sequences of filtered genes were constructed to validate the phylogenetic clades of nosZ against reference sequences. Sequences were aligned using MUSCLE (v3.8.1551)93 and trimmed with TrimAL (v1.4.1)94 with default settings. Maximum-likelihood trees were constructed with IQ-TREE (v2.2.0.3)95 with the “-m MFP -B 1000” options. The produced tree was visualized and beautified using the Interactive Tree of Life (iTOL; v6)96. Meanwhile, ESMFold97 was applied to predict the structure of each filtered gene. The 154 reference protein structures of NosZ were downloaded from AlphaFold DB98 and Protein Data Bank (PDB)99. A structural tree of NosZ was constructed using Foldtree (https://github.com/DessimozLab/fold_tree) based on a local structural alphabet100 and visualized using iTOL (v6)96. The structure of each gene in both the phylogenetic and structural trees was predicted using AlphaFold2101 and aligned against PDB using Foldseek (v8.ef4e960)102 with parameters “--tmscore-threshold 0.5 -e 0.001”. Cupredoxin-related protein active domains of nitrous-oxide reductase were investigated against the ECOD (Evolutionary Classification Of protein Domains) database103 with a TM-score >0.5 using the Foldseek easy-search module102.

Identification of nitric oxide dismutase gene (nod)
For the nod gene database search (see the workflow in Supplementary Fig. 4b), we initially downloaded reference protein sequences (n=1036) in NCBI’s databases to build an HMM model (available at https://doi.org/10.6084/m9.figshare.25650927). Then, nod genes in the non-redundant gene catalog and 3,164 MAGs were extracted using the above HMM model with hmmsearch in HMMER v3.3.2 using “-E 1e-5”. All potential nod genes shorter than 200 amino acids were excluded from further analysis. Following the methodology used for nosZ genes, these nod genes underwent verification through phylogenetic and structural trees. The reference protein structures for 44 Nod proteins were predicted using AlphaFold2101. The structures of the filtered nod genes were then predicted using AlphaFold2 and aligned with the reference protein structures using Foldseek (v8.ef4e960)102 with settings “--tmscore-threshold 0.5 -e 0.001”. In addition, diagnostic amino acid residues in the active center of the enzyme34 were identified in all proteins with a TM-score >0.5 using MAFFT (EMBL-EBI)104 and visualized with Jalview105.

Identification of hydrazine synthase and hydrazine dehydrogenase
To identify hzsA genes (see the workflow in Supplementary Fig. 4c), protein sequences in the non-redundant gene catalog and 3164 MAGs were firstly searched against NCycDB88, with the program DIAMOND blastp (version 2.0.14)89. Then, reference hzsA sequences (n=14) from the Greening lab metabolic marker gene databases90 were used to search for potential hzsA sequences in the non-redundant gene catalog and 3164 MAGs with DIAMOND blastp (version 2.0.14; --id 50)89. In addition, hzsA genes of non-redundant gene catalog and 3164 MAGs were extracted using HMMER v3.3.2 with the HMM profile “PF13486” from the InterPro database47. All identified hzsA genes were merged, and sequences shorter than 400 amino acids were filtered out. These genes were then verified by constructing phylogenetic and structural trees. The structure of each gene was predicted using AlphaFold2101 and aligned with PDB99 using Foldseek (v8.ef4e960)102. In addition, the active domain of the hydrazine synthase alpha subunit was searched in protein structures with a TM-score >0.5 using Foldseek’s easy-search module102.

To ensure the accuracy of identifying MAGs with hzsA genes, MAGpurify (v2.1.2) was used to detect contamination in MAGs through a combination of features and algorithms, including phylo-markers, clade-markers, tetra-freq, GC-content, and known-contam106. In addition, hzsB and hzsC sequences in MAGs were extracted with HMMER v3.3.2 using the “nitro.cycle.sub.hmm” model from the Metascan metabolic HMM database107. To further refine our search, the protein complexes of the previously undescribed clade hydrazine synthase (HzsABC) were predicted using AlphaFold (v2.0; model_preset=multimer)101. All structures were visualized and exported as images using PyMOL (http://www.pymol.org)108. In addition, the hzd gene associated with hydrazine dehydrogenase was identified in MAGs containing hzsA genes. To begin, we searched for reference protein sequences (n=49) in NCBI’s databases to construct the HMM model (available at https://doi.org/10.6084/m9.figshare.25650927). Subsequently, hzd genes were extracted using hmmsearch in HMMER, applying the “-E 1e-5” parameter. The structures were then predicted using AlphaFold2 and aligned with reference proteins in the PDB databases using Foldseek (v8.ef4e960) with settings “--tmscore-threshold 0.5 -e 0.001”.

MAG annotations and topological structure predictions
The MAGs were annotated using DRAM (v1.3.5)82 and Prokka (v1.14.6)109 with the default settings against KEGG, Pfam, MEROPS and dbCAN databases. Gene context was visualized using Chiplot (https://www.chiplot.online/), with the files produced by DRAM and Prokka as the input. DeepTMHMM (v1.0.24)110 and SignalP (v6.0)111 were employed to predict transmembrane topology and signal peptides of NosZ, Nod and HzsA proteins, respectively.

Identification of nifH, dsrA, mcrA and mobile genetic elements
To identify nifH genes associated with nitrogen fixation, protein sequences in the non-redundant gene catalog was first searched against NCycDB88, using the program DIAMOND (version 2.0.14)89 as mentioned above. Then, we utilized nifH reference sequences (n=1271) from the Greening lab metabolic marker gene databases90 to search for potential nifH sequences in the non-redundant gene catalog using DIAMOND blastp (version 2.0.14; --id 50)89. Subsequently, nifH genes were extracted using the NCBI’s Protein Family Models with HMM accession “TIGR01287.1” by HMMER v3.3.2. The genes were then verified through phylogenetic analysis as described in our previous study6. Conserved motifs (CXXR)6 were analyzed using MAFFT (EMBL-EBI)104 and visualized with Jalview105.

For the identification of dsrA and mcrA genes associated with sulfate reduction and anaerobic methane oxidation respectively, we employed eggNOG-mapper (v2.1.9; default parameters) for initial detection83,84. The reductive dsrA genes were further detected using DiSco (v1.0.0)112, while oxidative mcrA genes were confirmed via phylogenetic analysis as described in our previous study59. Classification of contigs belonging to mobile genetic elements (plasmids, proviruses, and viruses) was performed using Genomad v.1.5.0 with default parameters113.

Transcriptional activities of nitrogen loss genes
A total of 33 samples from various cold seeps—Haima, Qiongdongnan Basin, Shenhu area, and Jiaolong—were analyzed for metatranscriptomes, as outlined in our previous study114. Briefly, raw reads were quality filtered (parameters: --skip-bmtagger) using the Read_QC module within the metaWRAP (v1.3.2) pipeline115. SortMeRNA (v2.1)116 was employed to remove ribosomal RNA from quality-controlled reads using default settings. The abundances of transcripts for nosZ, nod, and hzsA genes were quantified by mapping filtered reads to a non-redundant gene catalog and genes of MAGs using Salmon (v.1.9.0; parameters: -validateMappings -meta)81. The results were expressed in transcripts per million reads (TPM)82.

Identification of nitrogen loss genes in surface sediments
A total of 34 samples from the Lingshui cold seeps in surface sediments, varying from 0 to 40 cmbsf, were analyzed in detail, as outlined in our previous study114. Briefly, the abundances of nosZ, nod and hzsA genes were quantified by mapping reads filtered using the Read_QC module within the metaWRAP to the non-redundant gene catalog and genes of 3164 MAGs using Salmon (v.1.9.0; parameters: -validateMappings -meta)81. The results were expressed in GPM.

Statistical analyses
Statistical analyses were performed using R v4.2.3. The normality of the data was assessed using Shapiro–Wilk tests prior to further analyses. The Kruskal-Wallis rank-sum test was used to compare the abundances of nitrogen loss genes across different depths, and their rates across different types of environments. Student’s t tests were utilized to examine differences in the relative abundance of nosZ, nod, and hzsA genes at varying depths, as well as nitrogen loss rates across non-seep and cold seep environments. Pairwise comparisons of environmental factors were conducted using Pearson’s correlation coefficients. In addition, relationships between the relative abundance of nosZ, nod, and hzsA genes and environmental factors was analyzed using Mantel tests.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "The authors bring an impressive suite of methods to bear on characterizing and evaluating the role of cold seep environments in the biologically mediated loss of fixed nitrogen from the ocean. They conclude that cold seeps are hot spots for N loss, mainly via denitrification. The significance of this hot spot phenomenon is measured in two ways: 1) Higher potential rates of N loss are detected in seep vs non-seep sites, based on data from three seep sites and one non-seep site and 2) Cold seeps are estimated to contribute 2% of the global N loss from marine sediments, while comprising only 0.087 % of the total ocean area.\nI don’t find either of these claims very compelling. 1) The data in figure 2 are convincing in showing that the newly measured rates at seep sites are significantly higher than the new measurements from non-seep sites. However, the number of non- seep measurements is very small. Secondly, as noted by the authors, the rates of both anammox and denitrification generally vary as a correlation of organic matter concentration or availability. That’s probably the main variable responsible for the difference observed here between seep and non-seep sites. I don’t think the scope of the new data is sufficient to support the global claim. 2) The measured rates are potential rates, which means they could be higher than in situ rates by a variable and unknown degree. Therefore it is really difficult to compare these rates to other reported rates and to know how any of them contribute to in situ rates. To compare the new rates to a global database of such rates would require assuming that the global data were distributed across a similar gradient of organic matter availability and it’s very unlikely that such a distribution can be validated for the available database of anammox and denitrification rates. Thirdly, the errors involved in extrapolating the rate data in the literature onto a global scale are surely larger than 2% (consider the uncertainty in both reported rates and the geographical calculations). Therefore, I don’t think we have the necessary precision to obtain a true estimate of the seep contribution.\nIn conclusion, this manuscript contains a lot of high quality data and solid findings. They don’t support the exaggerated claim of the title, however.\nSpecific comments\nL42: It’s hard to understand the significance of this number (6.16 Tg of N loss) without knowing what proportion of the total it represents. The structure of this sentence makes it sound like the number 6.16 refers the total N loss from seafloor surface sediments. Is this correct? We find out on L220 that this number is from seeps only, which is what I think you mean to intend here, but it’s not clear.\nL63: This sentence says that diazotrophs convert N2 into ammonium, nitrite and nitrate through biological N fixation. Strictly speaking, diazotrophs convert N2 into ammonium and other microbes convert the ammonium into nitrite and nitrate. Please reword.\nL80: The meaning of the sentence “Anammox, on the other hand….” Is not clear. “more energy efficient” than what? In what sense? The growth rates of anammox bacteria are usually pretty slow, much slower than that of conventional heterotrophic denitifiers, so that does not seem to be consistent with “more energy efficient”.\nWithout relying on oxygen? Neither of the two forms of denitrification described above “relies” on oxygen. One of them supposedly produces oxygen, but does not “rely” on the presence of oxygen.\nL90, L92: I don’t think it’s very useful to introduce the absolute rates of denitrification here. If the authors wish to compare rates they measure to previous reports from the literature, the place to do that would be where the rates from the present study are reported (in the results/discussion).\nL98: It’s not really appropriate to cite the Thamdrup et al. 2021 paper as evidence the anammox dominates in deep sea sediments. Thamdrup et al. did report strong dominance of anammox in the trench sediments, but the authors emphasized the anomalous nature of trench sediments compared to deep sea sediments. Those authors assumed that abyssal plain (non trench) sediments had lower rates of N loss overall due to the focusing of organic matter along the axis of the trench. It would be more accurate to replace “deep-sea” with “hadal” in this sentence, since hadal sediments do not represent deep- sea sediments and to say so misrepresents the meaning of the paper being cited.\nIn general, this paragraph contains information that is unnecessarily specific and overall not really necessary. If anammox and denitrification are important in marine sediments, they are likely occurring at cold seeps, no need to belabor the point.\nL162: The methods tell us that concentration of N2 gas was determined by GCMS. That doesn’t help me understand the units. 1.89 – 6.12 % of what? Total volume? Porewater volume? Above saturation? Please specify in the methods.\nL164: Why not help the reader understand the implications of the isotope data for N2 by saying: “The depleted values of d15N2 (-3.09 - 0.10 o/oo) indicated intense denitrification and anammox at Site F-14…. “ And that N2 concentration and d15N2 show approximately opposite depth distributions at this one station, which is what one would expect as evidence of denitrification.\nL168: Since the section is titled Results and Discussion, it would be appropriate to explain the significance of the TOC data, i.e., discuss it. Does higher TOC at the first three sites imply anything about difference in rates observed at those sites and Site F, and if so, what?\nFigure 2 is confusing. In panel a, data from many published sites are plotted. It looks like the rates of denitrification at the seep sites overlap with the rates in the deep sea, but the mean is slightly higher, whereas for the anammox, the mean is slightly lower. The only statistic appears to tell us that the rates differ among the five regions. Do the mean rates differ significantly between the seep sites and other deep-sea sites? That is the main comparison we need to make.\nThen in Figure 2b, seep and non-seep sites are compared. Do the non-seep sites include all the other sites in Figure 2a? I don’t think so. Therefore, it would be helpful to specify in the figure caption “(b) Potential N2 production rates measured in this study (denitrification and anammox) and the percent of N2 production attributed to anammox (ra, %) in non-seep and cold seeps sediments. P values for differences between non seep and cold seeps were computed using student's t-test.”\nOtherwise, it looks like seep sites have lower denitrification and anammox rates than most of the sediment environments shown in Figure 2a and does not make a strong case for the seep vs non seep difference.\nL208: The descriptive phrase “low temperature, high hydrostatic pressure, and stable carbon and energy supply” applies to most deep sea environments. I.e., it does not distinguish the adaptations needed to survive at seep vs non-seep environments.\nL 220: Given the uncertainty in rate measurements (and the fact that the new rates are potential rates) and the huge differences in area that a small change in assumptions could make, the 2% global contribution of N loss from seeps does not convince me that seeps are an important “hot spot”, either in the context of the global N budget or the uniqueness of seep environments.\nL 256 and elsewhere: Please specify either in the methods or at first use in the text that GPM (and TPM) means per million BASE PAIRS. One might logically think it meant genes per million genes.\nL 257: A difference in relative abundance of genes for N fixation and N loss does NOT translate to an imbalance in the flux of fixed nitrogen. Are there data for the rate of N fixation in deep sea or seep sediments?\nL 217: Does Fig. 9c refer only to the top 40cmbsf? This would imply that there are multiple depths within the 40 cmbsf dataset, but that is not mentioned or explained in the figure or text.\nL312: I’m not sure what definitive source of the known phylogenetic range for nosZ to which these data were compared, but at least some of the six new phyla listed in this line have been previously reported (Intrator et al 2024, Frontiers of Microbiology).\nFor the entire section from L279 – L429: This section contains a lot of interesting information, perhaps particularly the claim that the genetic capacity for anammox exists in several groups other than the Planctomycetes. However, it’s not always clear what the significance of the details of protein folding and motif binding have to do with the main story. I would suggest paring this section down to the significant essentials and moving the rest of it to another publication or at least the supplemental materials.\nL 494: I was surprised that the isotope tracer experiments were conducted with NO3 as the tracer (for determination of denitrification) and as the carrier/N pair for the anammox tracer 15NH4. Use of NO3 rather than NO2 could lead to an underestimate of both rates.\nL529: What completeness and contamination requirements were used to accept MAGs for further analysis?", "reviewer_2": "This paper explores the mechanisms of nitrogen loss in cold seep sediments, providing an important overview of this process, and suggesting global estimates of denitrification in seeps. The authors provide evidence for nitrogen loss, measure denitrification and anammox rates in slurry incubations, and investigate the taxonomic diversity of functional genes involved in these processes, using homology and structure-based approaches. This study is based on state-of-the-art techniques and difficult-to-collect samples, advancing our understanding of Earth’s biogeochemistry.\nThe manuscript can benefit from a better link between geochemistry and biochemistry. This can be achieved by better organization of the text, and deeper exploration of environmental correlations, as well as the taxonomic distribution of the studied genes – that is, answering the question of which taxa are the dominant denitrifiers, under which conditions?\nAdditional restructuring of the text is needed to clarify novelty and relevance, especially in parts that deal with the structures of proteins and gene clusters.\nSpecific comments:\nL51: Maybe remove “harsh conditions”? I would guess that the conditions elsewhere in the deep sea are also harsh, given, for example, the carbon limitation.\nL52-53: Please note that some seeps are not characterized by very low temperatures (e.g., the moderate temperatures in the eastern Mediterranean seeps, such as those from the Nile Deep Sea Fan).\nIntroduction: Please reorganize the introduction. In its current form, it moves back and forth from biogeochemical processes to enzyme biochemistry. For example, paragraph 2 introduces processes and respective enzymes, paragraph 3 deals with rates, and paragraph 4 concerns enzyme structure. There is no link between 3 and 4. This issue persists later in the manuscript, as it often struggles to link the exploration of enzymes\nResults and Discussion:\nSome paragraphs need to be adjusted to present the results coherently, in light of hypotheses and questions. For example, in lines 150-167, the paragraph may begin with: “Our geochemical measurements provide evidence for nitrogen loss in cold seep sediments.”. This pattern repeats throughout the text. Please try to revise the paragraph structure throughout the paper, so that the aim is explained early in each paragraph.\nL92 - “vary between 0 and 82.97…”\nL149 and elsewhere: please check if a hyphen is needed in “nitrogen loss”.\nL176 – Consider “Slurry experiments indicate substantial nitrogen loss rates in seep sediments”\nL189 – Consider “Although the anammox contributions were lower than those measured previously in deep-sea sediments…”\nL195 – The authors suggest that the positive correlation between denitrification and TOC suggests that heterotrophs are those who denitrify. Yet, the high TOC concentrations may be linked to higher methane fluxes, and thus enhanced chemosynthetic productivity. For example, some chemosynthetic Campylobacterota that often form seep biofilms appear to denitrify (e.g., Figure 3). Also, Supplementary Table 3 indicates that one of the most common NosZ-encoding organisms is the Beggiatoa. Moreover, substantial denitrification rates are often found near the sediment surface, where chemosynthesis is likely prominent.\nL226 – Maybe: “Near-surface sediments are hotspots of nitrogen loss gene diversity”\nLines 250+, 303+ and Fig.3 - It would be informative to see the per-taxa read abundance of the NosZ-encoding sequences, asking which lineages are the dominant denitrifiers.\nL245 – The MGE part is out of context. Please provide a working hypothesis to frame the MGE ideas.\nL255 – Please be careful when inferring process rates from gene read abundances. I suggest removing this part.\nL279 – The structural diversity is in genes/proteins, not in denitrifies – please rephrase\nL285-302 This part describes the general theory behind the discovery of NosZ proteins in the seep dataset. This is also the case of lines 324-336. It is not clear how this exploration is linked to the ecology/rates. Please either provide a better explanation or consider moving to the Supplement.\nL303 – I recommend exploring the diversity at lower taxonomic levels. Can you define the key denitrifying genera/species (as in the discussion of nod genes)? Can you assess their functionality (autotrophs/heterotrophs, etc..)?\nL324-336 – as mentioned above, it is difficult to trace the novelty in this paragraph. Are nos gene cluster structures in seep different from those described previously?\nL366 – “Genetic potential for ANNAMOX” is more suitable than “Annamox capabilities”.\nL402 – Discovery of ANNAMOX gene clusters in lineages other than those described previously requires careful examination, beyond the completeness stats of the bins. To make the statement stronger, please check for uniform coverage of contigs. Also, please examine the taxonomic distribution (e.g., using BLAST against NCBI’s nr) of the genes described, as well as the adjacent ones (such as those that encode the upstream genes related to c-type 406 cytochromes and Ca-activated chloride channels).\nL409 – What do you mean by “traditional” and “novel” HzsA? Please define earlier in the text.\nL542 – Markov\nSupplementary Table 8 – How were the abundance percentages calculated? They seem very low…\nPlease add legends to the Supplementary Tables in the main supplementary document (the one with the figures). For example, it is not clear which units are used in Tables S4 and S5.\nSupplementary Table 24 – Please provide the stats for all the MAGs used in this study."}, "author_response": "Response: We sincerely thank the reviewer for the careful evaluation of our work and for the insightful comments provided. We have carefully considered all the points raised and have responded to each issue in detail below.\n[...]\nResponse: We acknowledge that the initial number of non-seep measurements in our study was limited. To address this, we have expanded our dataset to include nitrogen loss rates from additional non-seep sediment samples (five additional cores, totaling 23 samples) from the Qiongdong, Shenhu, Lingshui and Site F regions. Furthermore, we supplemented our dataset with nitrogen loss rates from cold seep sediment samples (four additional cores, totaling 23 samples) from the Lingshui, Haima and Site F regions. These additions provide a broader spatial representation and strengthen the comparative analysis between seep and non-seep sites.\n[...]\nResponse: We agree that organic matter concentration is likely a key factor driving the observed differences in nitrogen loss rates between seep and non-seep sites. Our results showed that TOC levels at some cold seep sites (Lingshui-6/-10/-11/-12, Haima-7/-8/-9, and Site F-14/-15) were much higher than other sites, suggesting that these sites could provide more electron donors for heterotrophic denitrification and leading to potentially higher nitrogen loss rates (Supplementary. Fig. 2). This finding supports TOC as a primary factor influencing denitrification in cold seep sediments. Additionally, as Reviewer 2 noted, chemosynthetic microorganisms enriched in cold seep sediments likely play an important role in the denitrification process. For instance, certain chemosynthetic Campylobacterota, known to form biofilms in seep environments, are capable of denitrification as shown in Figure 3. Other sulfur-oxidizing bacteria, such as Beggiatoaceae, which is abundant in cold seep sediments and encodes the nosZ gene (Supplementary Table 3), may also contribute significantly to nitrogen cycling in these habitats (L194–L200).\n[...]\nResponse: We recognize the limitations of using potential rather than in situ rates, as slurry incubations may yield nitrogen loss rates that exceed actual in situ values. In the revised manuscript, we have clarified that our rates represent potential values. Additionally, we acknowledge that extrapolating these findings on a global scale involves significant uncertainty, especially given the variability in cold seep sizes and distributions. In response, we conducted a more detailed literature review on global cold seep areas and found that the estimated global cold seep area likely ranges between 2.05 and 3.15 × 105 km². Based on this updated area, we recalculated the global nitrogen loss flux from denitrification and anammox. We also agree that comparing our results to a global database assumes a similar gradient of organic matter availability, which is difficult to validate. To address this, we have added a cautionary note in the manuscript, clarifying that our global estimates are approximations and that variations in organic matter concentrations across regions limit direct comparisons.\n[...]\nResponse: We agree that the original title may have overstated the findings. In response, we have revised the title to better align with the data and the study’s exploratory nature, addressing the limitations related to e.g. potential rather than in situ nitrogen loss rates and the extrapolation challenges. By describing cold seeps as potential hotspots, we emphasize their likely significance without overstating the findings.\nThe revised title now reads: “Cold seeps are potential hotspots of deep-sea nitrogen loss driven by microorganisms across 21 phyla”.\n[...]\nResponse: We thank the reviewer for detailed reading and insightful comments. Detailed responses for each comment are provided below.\nThe manuscript can benefit from a better link between geochemistry and biochemistry. This can be achieved by better organization of the text, and deeper exploration of environmental correlations, as well as the taxonomic distribution of the studied genes – that is, answering the question of which taxa are the dominant denitrifiers, under which conditions?\nResponse: We appreciate your insightful comment and have revised the manuscript to directly address the link between geochemistry and biochemistry. To achieve this, we reorganized the text to better integrate environmental correlations with the taxonomic distribution of nitrogen loss genes, explicitly identifying the dominant denitrifiers and the conditions under which they thrive in cold seep environments."}, {"round": 2, "reviewer_comments": {"reviewer_1": "The authors have thoroughly and meticulously revised the manuscript in response to my suggestions and queries. Importantly, they have provided substantial additional data and results, which do much to strengthen the conclusions and improve the robustness of their statistical comparisons. They have also toned down the claims about the global significance of seeps as a site for loss of fixed nitrogen. My main concern in making that claim is the large amount of error/uncertainty involved in scaling up the rate measurements and area estimate. These uncertainties mean that a contribution on the order of 2% is well within the uncertainty of the estimate of the total value, so it’s not reasonable to say that seeps make an important contribution to the global N loss budget.\nHowever, that does not remove the global implications entirely. This expanded investigation of global seeps does show convincingly that seeps are something of a “hot spot” in comparison to the rest of the deep sea. So, on a global basis, that’s a robust finding.\nThe additional samples have been incorporated into the analysis, and in the revised figures. However, it looks like the revised Figure 1 did not make it into the manuscript file. The version of Figure 1 in the rebuttal document contains new stations (a total of 45 sites are listed in the legend) but the version of the figure in the actual manuscript looks identical to the original and contains only 36 sites in the legend.\nLine 585: “Details” not “Detailed”", "reviewer_2": "The authors adequately addressed my queries. I recommend accepting the manuscript."}, "author_response": "Response: We are very pleased to hear this, and we appreciate the reviewer’s helpful comments and corrections.\n[...]\nResponse: We thank the reviewer for reading our paper carefully. We have now corrected version of Figure 1 which includes the additional stations (45 sites listed in the legend) in the revised manuscript.\n[...]\nResponse: Thank you for your valuable feedback on this issue. We think this comment pertains to “Detailed”, mentioned in L331 of the original manuscript. Corrected as suggested.\n[...]\nResponse: We are very pleased to hear this, and we thank the reviewer for constructive comments on the manuscript."}]

Decision letter