Diversity and N2O Production Potential of Fungi in an Oceanic Oxygen Minimum Zone

Fungi in terrestrial environments are known to play a key role in carbon and nitrogen biogeochemistry and exhibit high diversity. In contrast, the diversity and function of fungi in the ocean has remained underexplored and largely neglected. In the eastern tropical North Pacific oxygen minimum zone, we examined the fungal diversity by sequencing the internal transcribed spacer region 2 (ITS2) and mining a metagenome dataset collected from the same region. Additionally, we coupled 15N-tracer experiments with a selective inhibition method to determine the potential contribution of marine fungi to nitrous oxide (N2O) production. Fungal communities evaluated by ITS2 sequencing were dominated by the phyla Basidiomycota and Ascomycota at most depths. However, the metagenome dataset showed that about one third of the fungal community belong to early-diverging phyla. Fungal N2O production rates peaked at the oxic–anoxic interface of the water column, and when integrated from the oxycline to the top of the anoxic depths, fungi accounted for 18–22% of total N2O production. Our findings highlight the limitation of ITS-based methods typically used to investigate terrestrial fungal diversity and indicate that fungi may play an active role in marine nitrogen cycling.


Introduction
Oceanic oxygen minimum zones (OMZs) are characterized by a sharp oxycline and redox gradient in the water column [1]. As a result, OMZs support diverse microbial communities that directly impact the global biogeochemical cycling of nitrogen, carbon, sulfur, and trace metals [2][3][4][5][6]. As in many other types of marine environments, bacteria and archaea have been the focus of microbial ecology research in OMZs, whereas microbial eukaryotes, in particular fungi, have received much less attention [7,8].
An early cultivation-based survey of marine fungi in the Indian Ocean that included the Arabian Sea OMZ found Rhodotorula rubra and Candida atmosphaerica to be cosmopolitan, and the yeast population densities ranged from 0-513 cells per liter of seawater [9]. In the eastern tropical South Pacific (ETSP) OMZ off the coast of Chile, high summertime fungal biomass in the water column have been reported [10], including diatom parasites from the phylum Chytridiomycota [11]. In the third and the largest open ocean OMZ, the eastern tropical North Pacific (ETNP), protists diversity has been investigated by sequencing the V4 region of 18S small subunit rRNA genes [7], but little is known about the diversity and function of fungi in this environment.
Fungi in the water column are generally thought to contribute to organic matter recycling, particularly in particle-associated environments [12][13][14]. High hydrolytic activity on proteinaceous substrates in large size fractions (>25 µm and >90 µm) have been reported in the water column of the ETSP and attributed to fungi given low bacterial biomass in those size fractions [15]. However, it remains unclear if fungi in the water column of OMZs play a role in nitrogen cycling. Discovered in the early 1990s, fungal denitrification is known as a process that reduces nitrate or nitrite with nitrous oxide (N 2 O) as the end-product [16,17]. This adds to the multiple other pathways and processes (ammonia oxidation, bacterial denitrification, and chemodenitrification) that can produce N 2 O [18], a potent greenhouse gas and ozone-depleting agent [19]. Many fungal strains have been found to have the ability to produce N 2 O [20], including an Aspergillus terreus strain isolated from the Arabian Sea OMZ [21]. In marine environments, fungal denitrification with N 2 O as the end-product has been reported from coastal marine sediment in India and Germany [22,23], but its potential contribution in the water column remains unclear.
We investigated the fungal community composition in the eastern tropical North Pacific oxygen minimum zone by sequencing the internal transcribed spacer region 2 (ITS2) and classifying shotgun metagenome reads. To estimate the fungal contribution to N 2 O production in the water column, we used a selective inhibition method combined with 15 N-labeled tracer incubation experiments. Fungal communities evaluated by ITS2 sequencing were dominated by the phyla Basidiomycota and Ascomycota at most depths. The metagenome dataset showed that early-diverging fungi accounted for about one third of the fungal community, and the subsurface peaks of fungal abundance coincided with both cyanobacterial abundance and eukaryotic algal abundance. Incubation experiments suggest a possible role of fungi in N 2 O production in the water column of the ETNP OMZ.

Site Description and Seawater Filtration
In March 2018, aboard the R/V Sally Ride in the eastern tropical North Pacific oxygen minimum zone, two stations were visited to study fungal diversity and the potential fungal contribution to N 2 O production ( Figure 1a). Dissolved oxygen concentration was determined using the SBE 43 dissolved oxygen sensor attached to the conductivity, temperature, and depth (CTD) rosette. Seawater was collected at multiple depths spanning from the oxycline to the anoxic depths ( Figure 1b)  biomass in those size fractions [15]. However, it remains unclear if fungi in the water column of OMZs play a role in nitrogen cycling. Discovered in the early 1990s, fungal denitrification is known as a process that reduces nitrate or nitrite with nitrous oxide (N2O) as the end-product [16,17]. This adds to the multiple other pathways and processes (ammonia oxidation, bacterial denitrification, and chemodenitrification) that can produce N2O [18], a potent greenhouse gas and ozone-depleting agent [19]. Many fungal strains have been found to have the ability to produce N2O [20], including an Aspergillus terreus strain isolated from the Arabian Sea OMZ [21]. In marine environments, fungal denitrification with N2O as the end-product has been reported from coastal marine sediment in India and Germany [22,23], but its potential contribution in the water column remains unclear. We investigated the fungal community composition in the eastern tropical North Pacific oxygen minimum zone by sequencing the internal transcribed spacer region 2 (ITS2) and classifying shotgun metagenome reads. To estimate the fungal contribution to N2O production in the water column, we used a selective inhibition method combined with 15 N-labeled tracer incubation experiments. Fungal communities evaluated by ITS2 sequencing were dominated by the phyla Basidiomycota and Ascomycota at most depths. The metagenome dataset showed that early-diverging fungi accounted for about one third of the fungal community, and the subsurface peaks of fungal abundance coincided with both cyanobacterial abundance and eukaryotic algal abundance. Incubation experiments suggest a possible role of fungi in N2O production in the water column of the ETNP OMZ.

Site Description and Seawater Filtration
In March 2018, aboard the R/V Sally Ride in the eastern tropical North Pacific oxygen minimum zone, two stations were visited to study fungal diversity and the potential fungal contribution to N2O production ( Figure 1a). Dissolved oxygen concentration was determined using the SBE 43 dissolved oxygen sensor attached to the conductivity, temperature, and depth (CTD) rosette. Seawater was collected at multiple depths spanning from the oxycline to the anoxic depths ( Figure 1b [24]. (b) Depths sampled for nitrous oxide (N2O) production experiments are marked by filled symbols. Color contour shows oxygen concentrations measured during this cruise using a Seabird SBE 43 dissolved oxygen sensor.
To collect particulate material at different size fractions, seawater was sequentially filtered through a 47 mm Whatman Grade 541 acid-hardened cellulose filter paper (22 μm nominal particle retention rating, GE Healthcare 1541-047, Marlborough, MA, USA), a 47 mm polycarbonate filter (2.0 μm nominal pore size, Millipore Isopore TTTP-04700, To collect particulate material at different size fractions, seawater was sequentially filtered through a 47 mm Whatman Grade 541 acid-hardened cellulose filter paper (22 µm nominal particle retention rating, GE Healthcare 1541-047, Marlborough, MA, USA), a 47 mm polycarbonate filter (2.0 µm nominal pore size, Millipore Isopore TTTP-04700, Burlington, MA, USA), and a Sterivex filter (0.22 µm nominal pore size, Millipore SVGP01050, Burlington, MA, USA), using a peristaltic pump filtration at a flow rate < 50 mL/min. For each sample set, 23 to 55 L of seawater was filtered (Table S1). Each 47 mm filter was stored in a 47 mm petri dish and flash-frozen in liquid nitrogen before storage at −80 • C.

DNA Extraction
In the laboratory, DNA was extracted using the DNeasy Plant Mini Kit (QIAGEN Cat No. 69104, Germantown, MD, USA) following the DNeasy Plant Handbook [25], except that the cell disruption step was customized for our samples. Filter paper from the Sterivex filters were extracted from the plastic case using a snap-blade knife. All filters were first cut into 2 by 2 mm pieces using sterilized scissors and transferred into 2 mL screw cap tubes containing 1 mL of 0.5 mm zirconia/silica beads (Biospec products #11079105z, Bartlesville, OK, USA), 600 mL of buffer AP1, and 6 µL of RNase A. Bead beating of the samples was performed for 90 s using a Biospec Mini-BeadBeater-16 and was followed by incubation at 65 • C for 10 min. After centrifugation at 20,000× g for 5 min, the supernatant in each sample tube was transferred to a fresh 2 mL microcentrifuge tube and neutralized with 195 µL of Buffer P3. The remaining DNA extraction steps followed the DNeasy Plant Handbook without modifications. An extraction blank was included for each batch of extraction procedure. DNA yield was quantified using a Qubit fluorometer (ThermoFisher Scientific, Waltham, MA, USA) following the manufacturer's instructions. All extracted DNA samples were stored at −80 • C until amplicon library construction.

Sequencing and Analysis of the ITS2 Region
The second region of the internal transcribed spacer (ITS2) with flanking regions in the 5.8S and 28S ribosomal RNA was targeted for amplicon library preparation following the Illumina 16S Metagenomic Sequencing Library Preparation [26] with the following modifications. The amplicon PCR was performed using primers ITS3tagmix [27] and ITS4tag001 [28] (Table S2)  Raw sequence reads were first trimmed using ITSxpress to contain only the ITS2 region [29]. Trimmed reads were merged, quality filtered, dereplicated, and denoised following the USEARCH pipeline [30] to generate amplicon sequence variants (ASVs). Taxonomic assignment was performed using a combination of Naïve Bayes classifier implemented in QIIME2 [31] and BLASTn [32] against full UNITE + INSD dataset v02.02.2019 [33] and curated manually. To determine putative fungal denitrifiers in our samples, we performed a closed reference OTU picking [34] of our ASVs and ITS2 sequences from fungi tested for N 2 O production [20] using UCLUST [35] implemented in QIIME v1.9.1 [36] against the UNITE database [33]. Raw reads generated in this study are available at the National Center for Biotechnology Information (NCBI) under BioProject PRJNA623945.

Analysis of Fungal Diversity and Function from Metagenomes
As an independent approach to evaluate the fungal diversity in the eastern tropical North Pacific oxygen minimum zone, we investigated a metagenome dataset sampled at a nearby station in March 2012 when the hydrographic conditions were highly similar ( Figure S1) [37]. Raw reads were first filtered using the tool BBDuk Version 38.73 [38] with the options "ktrim=r ordered minlen=51 minlenfraction=0.33 mink=11 tbo tpe rcomp=f k=23 ftm=5". Adapters were trimmed from the BBDuk-filtered reads using the tool Trimmomatic Version 0.39 [39] with the options "ILLUMINACLIP:$adapters:2:30:10 LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 MINLEN:100". Reads that passed both quality filtering and adapter trimming were queried against the NCBI nr database using DIA-MOND [40] with an e-value threshold of 1 × 10 −5 and the options "-sensitive -min-orf 20". The resultant NCBI taxonomy ID was used to assign taxonomy to each read. In search for the presence of the gene diagnostic for fungal N 2 O production, the cytochrome P450 nitric oxide reductase (P450nor) [41], we searched the metagenome assemblies under the same BioProject (PRJNA350692) for putative P450nor using HMMER v3.2.1 [42] against a P450nor profile, which includes both eukaryotic and prokaryotic cytochrome P450 genes [43]. All putative P450nor hits were checked against the NCBI nr database using blastp [32] to determine taxonomy.

Measurements of Potential N 2 O Production Rates
To measure potential rates of N 2 O production, parallel incubation experiments were performed by adding 0.1 mL of 5 mM 99% pure 15 N-labeled potassium nitrate or 0.1 mL of 0.8 mM 99% pure 15 N-labeled ammonium chloride (Cambridge Isotopes, Cambridge, MA, USA) to 120 mL glass serum bottles containing freshly collected seawater. To minimize the introduction of atmospheric oxygen, each bottle was overflown three times its volume with water directly from Niskin bottles before filling and crimp-sealing it with grey butyl rubber stopper and aluminum caps. Headspace was created in each bottle by replacing 2 mL of seawater with ultra-high pure helium (Airgas HE UHP300, Radnor, PA, USA). End points were taken by adding 1 mL of 50% (w/v) zinc chloride to separate parallel incubations approximately 0 and 24 h after incubations began at Station 1 and 0, 12, and 24 h after incubations began at Station 2. The concentration of ammonium in seawater was measured onboard according to the fluorometric method of Holmes et al. [44], with a detection limit of 15 nmol L −1 . Nitrate concentrations were assayed using a Lachat Flow Injection Analyzer at the Analytical lab at the Marine Science Institute, University of California, Santa Barbara following standard analytical methods [45]. The detection limit of nitrate was 0.2 µmol L −1 .
The potential contribution of fungal N 2 O production to total N 2 O production was determined by incubations with 15 N-labeled nitrate ( 15 NO 3 − ) and chloramphenicol (87.7 mg L −1 final concentration), which was applied to each incubation one hour before the addition of 15 N tracers to inhibit all prokaryotic activities. All solutions added to the incubation bottles were purged by ultra-high pure helium for one hour at a flow rate of 40 mL min −1 . The difference between fungal N 2 O production and total N 2 O production from incubations with 15 NO 3 − is attributed to bacterial denitrification ( Figure S2). N 2 O production rates measured in incubations with 15 N-labeled ammonium ( 15 NH 4 + ) are attributed to archaea and/or bacterial nitrification.
The quantity and isotopic composition of dissolved N 2 O was determined using a Delta XP isotope ratio mass spectrometer coupled to a purge-and-trap front end. The detection limit was 1.0 nmol N, and the precision for δ 15 N was 2.0‰ (n ≥ 3). The rate of N 2 O production (R N2O ) was calculated from the equation [46]:  Figure S4).

Fungal Diversity Assessed by Sequencing the ITS2 Region
Assessment of the fungal community in the eastern tropical North Pacific using the internal transcribed spacer region 2 (ITS2) revealed that taxa from the phyla Basidiomycota and Ascomycota dominated at most depths at both stations ( Figure 2). The relative abundance of Basidiomycota was higher than Ascomycota at nearly all depths (Table S3). The most prevalent and abundant taxon is the basidiomycetous yeast family Sporidiobolaceae, primarily consisting of the genera Rhodotorula, Rhodosporidiobolus, and Sporobolomyces (Table S4). Sporidiobolaceae tend to have a higher relative abundance in the 0.2-2 µm size fraction. On the other hand, Aureobasidium (Ascomycota) and Exobasidiomycetes (Basidiomycota, primarily Meira) were enriched in the larger size fractions (2-22 and >22 µm). In contrast, the basidiomycetous yeast Malassezia, when detected, were enriched only in the 2-22 µm size fraction. At Station 2, most of the fungal community cannot be classified even at the phylum level based on ITS2 sequences, indicating the presence of novel fungal lineages in the oxycline of ETNP oxygen minimum zone.

Fungal Diversity Assessed from Metagenomes
Fungal community composition assessed by metagenomic reads also showed the dominance of Dikarya fungi, but in contrast to the results from ITS2 sequencing, the relative abundance of Ascomycota was consistently higher than that of Basidiomycota ( Figure  3). Surprisingly, over one third of the fungal community belong to early-diverging phyla

Fungal Diversity Assessed from Metagenomes
Fungal community composition assessed by metagenomic reads also showed the dominance of Dikarya fungi, but in contrast to the results from ITS2 sequencing, the relative abundance of Ascomycota was consistently higher than that of Basidiomycota ( Figure 3). Surprisingly, over one third of the fungal community belong to early-diverging phyla including Mucoromycota, Zoopagomycota, Chytridiomycota, Blastocladiomycota, Cryptomycota, and Microsporidia. The fungal community composition from 60 to 300 m was uniform, while there was a trend of increasing relative abundance for Ascomycota with depth.

Relative Abundance of Fungi Compared to Other Taxa
Mining the metagenomes allowed us to estimate the relative abundance of fungi and other taxa as part of the overall microbial community. In each metagenome, 23-57% of the reads had a positive hit against the NCBI nr database with an e-value threshold of 1 × 10 −5 ( Figure S5a). Fungi accounted for 0.02-0.22% of all classifiable reads, with a subsurface peak at 70 m (Figure 4a). The relative abundance of fungal reads decreased with depth below 70 m but showed a small increase at 140 m, which was the top of the oxygen deficient layer of the water column. Ciliophora, Oomycetes, and Dinophyceae were similar to fungi in both abundance and distribution. The abundance of reads classified as eukaryotic algae also showed a subsurface peak at 70 m (Figure 4b). The abundance of cyanobacteria decreased with depth overall, but there was an increase at 110 m, coinciding with the deep chlorophyll maximum situated immediately above the anoxic depths.

Relative Abundance of Fungi Compared to Other Taxa
Mining the metagenomes allowed us to estimate the relative abundance of fungi and other taxa as part of the overall microbial community. In each metagenome, 23-57% of the reads had a positive hit against the NCBI nr database with an e-value threshold of 1 × 10 −5 ( Figure S5a). Fungi accounted for 0.02-0.22% of all classifiable reads, with a subsurface peak at 70 m (Figure 4a). The relative abundance of fungal reads decreased with depth below 70 m but showed a small increase at 140 m, which was the top of the oxygen deficient layer of the water column. Ciliophora, Oomycetes, and Dinophyceae were similar to fungi in both abundance and distribution. The abundance of reads classified as eukaryotic algae also showed a subsurface peak at 70 m (Figure 4b). The abundance of cyanobacteria decreased with depth overall, but there was an increase at 110 m, coinciding with the deep chlorophyll maximum situated immediately above the anoxic depths.

Potential Contribution of Fungal N2O Production
N2O production rates from incubation experiments with either 15 NO3 − or 15 NH4 + were detected at multiple oxycline depths of the ETNP oxygen minimum zone (OMZ), and the maximum rate was found at the oxic-anoxic interface ( Figure 5). The rates of N2O production from incubations with 15 NO3 − and chloramphenicol were used as an approximation of fungal N2O production, and they ranged from 0% of total N2O production from 15 NO3 − at 275 m at Station 2 to 56% of total N2O production from 15 NO3 − at 90 m at Station 1. N2O production from incubations with 15 NO3 − , both with and without chloramphenicol, was lower at elevated in situ oxygen (O2) concentration. When integrated from the oxycline (60 m at Station 1 and 90 m at Station 2) to the oxic-anoxic interface, fungal N2O production approximated by incubations with 15 NO3 − , and chloramphenicol accounted for 18-22% of total N2O production (Table S5).

Potential Contribution of Fungal N 2 O Production
N 2 O production rates from incubation experiments with either 15 NO 3 − or 15 NH 4 + were detected at multiple oxycline depths of the ETNP oxygen minimum zone (OMZ), and the maximum rate was found at the oxic-anoxic interface ( Figure 5). The rates of N 2 O production from incubations with 15 NO 3 − and chloramphenicol were used as an approximation of fungal N 2 O production, and they ranged from 0% of total N 2 O production from 15 (Table S5).

Search for Fungal Denitrifiers and Functional Genes
To search for putative fungal denitrifiers in the water column of the eastern tropical North Pacific OMZ, we first analyzed the ITS2 amplicon dataset generated in this study. Of all 237 ASVs, only one (ASV211) was clustered at 97% similarity level with the ITS2 sequences from fungal strains previously shown to produce N2O (Chaetomium sp.) [20]; this ASV was present only in the 0.2-2 μm size fraction at 83 m from Station 1 and at a relative abundance of 0.2% (Table S6). However, one ASV (ASV8) was 96.95% similar to the ITS2 region of the N2O-producing Penicillium melinii [20], which was isolated from seawater [48]. The relative abundance of ASV8 ranged from 0.1 to 8.5% at hypoxic and

Search for Fungal Denitrifiers and Functional Genes
To search for putative fungal denitrifiers in the water column of the eastern tropical North Pacific OMZ, we first analyzed the ITS2 amplicon dataset generated in this study. Of all 237 ASVs, only one (ASV211) was clustered at 97% similarity level with the ITS2 sequences from fungal strains previously shown to produce N 2 O (Chaetomium sp.) [20]; this ASV was present only in the 0.2-2 µm size fraction at 83 m from Station 1 and at a relative abundance of 0.2% (Table S6). However, one ASV (ASV8) was 96.95% similar to the ITS2 region of the N 2 O-producing Penicillium melinii [20], which was isolated from seawater [48]. The relative abundance of ASV8 ranged from 0.1 to 8.5% at hypoxic and anoxic depths (Table S6), and it belongs to the same family (Aspergillaceae) as the N 2 Oproducing fungal strain isolated from the Arabian Sea OMZ [21]. While this putative denitrifier could potentially produce N 2 O in the water column of the OMZ, its potential contribution to the total N 2 O production is estimated to be less than 0.3% given the low N 2 O yield by P. melinii in the laboratory (Appendix A). Continuing the search for putative fungal denitrifiers, we analyzed a previously published March 2012 metagenome dataset from a nearby station [37] for the occurrence of P450nor. All P450nor hits identified by a hidden Markov model (HMM) profile search [43] were prokaryotic genes, most of which from the phylum Actinobacteria (Table S7).

Fungal Diversity in Different Size Fractions
In the eastern Tropical North Pacific (ETNP) oxygen minimum zone, the most prevalent and abundant fungal taxa in our ITS2-based survey was the basidiomycetous yeast from the family Sporidiobolaceae, which are usually characterized by the production of carotenoids that color their cells red, pink, or orange and, hence, have the name "red yeasts". The high relative abundance of red yeasts in the 0.2-2 µm size fraction ( Figure 2) is consistent with their typically small cell size [49]. A previous study on yeasts in the Indian Ocean, which includes the Arabian Sea oxygen minimum zone in the northern part, has also found red yeasts to be the predominant taxa [9]. On the other hand, Aureobasidium (Ascomycota) and Exobasidiomycetes (Basidiomycota, primarily Meira) were enriched in the larger size fractions (2-22 µm and >22 µm), consistent with a larger size [50,51], and indicating a likely association with particles. In contrast, the basidiomycetous yeast Malassezia, when detected, were enriched only in the 2-22 µm size fraction, suggesting that they were not associated with particles larger than 22 µm.
There is increasing recognition of marine fungi as a key component of the marine microbiome and biogeochemical cycles [12,13,15,52,53], but our knowledge about their diversity and function is far less compared to other microbial eukaryotes [54]. This is particularly the case in the water column of the open ocean, where metabarcoding surveys of fungal diversity are often unable to classify most of the fungal community [55,56]. In one of the two stations sampled in this study, we also could not resolve the taxonomy for most of the fungal community (Figure 2), highlighting the limitation of an ITS-based approach to study fungal diversity. Because the public databases for ITS sequences are primarily based on studies of terrestrial environments or fungal strains, this suggests that the open ocean water column harbors previously undiscovered lineages of fungi.

Ecology of Marine Fungi in the Oxygen Minimum Zone
By classifying metagenome reads, we showed that about a third of the fungal community were from early-diverging lineages [57], including Mucoromycota, Zoopagomycota, Chytridiomycota, Blastocladiomycota, Cryptomycota, and Microsporidia ( Figure 3). Such diversity was previously undiscovered due to both the difficulty in cultivation and the consequent lack of representation in ITS databases. It is unclear what ecological roles these early-diverging fungi play, except for Chytridiomycota, which are typically associated with phytoplankton such as diatoms [58]. Nevertheless, the depth profile of the relative abundance of fungi in relation to other microbial taxa provides a hint (Figure 4). The cooccurrence of subsurface peaks (at 70 m) of the relative abundance of fungi and eukaryotic algae suggests that fungi are directly associated with eukaryotic algae, perhaps as parasites. A secondary peak of the relative abundance of fungi (at 140 m) was observed below that of the deep chlorophyll maximum (at 110 m) typically found at the oxic-anoxic interface of oxygen deficient waters [59]. Hence, the fungal communities at the secondary peak (at 140 m) are likely predominantly saprotrophic, feeding on the particles formed at the base of the mixed layer.
Surprisingly, we did not observe a pronounced effect of oxygen concentration on the fungal community composition, evaluated by either ITS2 amplicon sequencing or metagenomes. In contrast, the protist community in the anoxic depth of the ETNP and eastern tropical South Pacific oxygen minimum zone was enriched in Syndiniales, euglenozoan flagellates, and acantharean radiolarians [7,8]. This may be a result of the versatile respiratory/fermentative metabolisms fungi possess, but it could also be due to the inability of the methods used in this study to detect novel fungal lineages from the oxygen deficient waters. The metagenome-based approach avoids the typical biases associated with amplicon sequencing such as primer bias, but it is limited by the sequences available in the chosen database (NCBI nt in this study).

Fungal N 2 O Production in the Oxygen Minimum Zone
In order to distinguish the N 2 O production by fungi from bacteria and archaea, we combined 15 N tracer incubation experiments with selective inhibition using the antibiotic chloramphenicol ( Figure S2). Chloramphenicol is a broad-spectrum antibiotic that has been used to isolate anaerobic gut fungi from the rumen microbiome [60], and its final concentration used in our incubations was scaled by the typical bacterial cell density in the rumen vs. seawater. Nonetheless, we did not make direct measurements of the specificity and effectiveness of chloramphenicol on inhibiting bacterial and archaeal activities in seawater. Therefore, it is possible that the N 2 O production rates we measured from incubations with 15 NO 3 − and chloramphenicol include N 2 O from partially inhibited bacterial denitrification. Consequently, we conservatively interpret those rates as an upper limit of potential fungal N 2 O production.
When integrated from the oxycline to the oxic-anoxic interface, fungal N 2 O production accounted for 18-22% of total N 2 O production (Table S5). While we interpret this as the maximum possible contribution of fungi to N 2 O production, we suggest that fungal denitrification could be an important pathway for N 2 O production in oceanic oxygen minimum zones. Denitrification, the sequential reduction in NO 3 − to N 2 , is known to be inhibited by trace amounts of O 2 [61]. In OMZs, it was demonstrated that 297 nM of O 2 repressed 50% of total N 2 O production at the oxic-anoxic interface [62]. In this study, N 2 O produced in incubation with 15 NO 3 − (primarily via denitrification) was inhibited by increasing levels of in situ O 2 concentration. However, the inhibitory effect of O 2 appeared to be less pronounced on fungal denitrification than on bacterial denitrification at low levels, revealing a potential niche (0.0 < O 2 < 0.93 µM) for fungi capable of denitrification. Since this potential niche is shallower than the oxygen deficient waters, the potential fungal N 2 O production can have a higher chance of reaching the ocean-atmosphere interface than N 2 O produced at deeper depths. Alternatively, N 2 O produced by fungi can be taken up by bacteria with the atypical ("Clade II") nitrous oxide reductase [63] (primarily Flavobacteria and Chloroflexi), which was shown to have a peak in relative abundance at depths immediately above the oxic-anoxic interface [37].

Molecular Evidence for Fungal N 2 O Production
The search for putative fungal denitrifiers using ITS2 sequence identity suggests an extremely low abundance of known N 2 O-producing fungi in the water column of oxygen minimum zones. However, it should be noted that the collection of N 2 O-producing fungi used to identify these putative denitrifiers consists of soil fungi exclusively [20], so fungal lineages capable of N 2 O production from the open ocean were likely excluded. Therefore, there may be other N 2 O-producing fungi from the ETNP OMZ unidentified by this approach.
The absence of fungal P450nor in the metagenomic dataset we queried may be attributed to insufficient sequencing depth, given the low percentage (0.02-0.22%) of fungal reads classified by DIAMOND against the NCBI nr database [40] (Figure 4a). Even under the most simplistic assumption that there was only one fungal species present and it possessed P450nor in its genome, the sequencing depth in most samples was insufficient to recover just one copy of fungal P450nor (Table S8). Additionally, most DNA extraction protocols applied in published metagenome studies (including [37]) are not customized for disrupting chitinous cell walls. This likely resulted in under-sampling DNA from fungi, of which the biomass is low in the ocean water column (0.01-0.12 µg of carbon L −1 ) [10], especially compared to bacteria (estimated average of 10.5 µg of carbon L −1 ) [64]. Finally, it should be noted that the detection of fungal P450nor genes in metagenomes does not necessarily imply fungal denitrification, as P450nor in certain fungal genomes appear to be involved in secondary metabolisms instead [43].

Conclusions
Our findings highlight the previously unrecognized fungal diversity in the eastern tropical North Pacific oxygen minimum zone, particularly from the early-diverging taxa as revealed by analysis of shotgun metagenomes. The depth distribution pattern of fungi in relation to cyanobacteria and eukaryotic algae suggest direct association in the mixed layer of the water column and indirect feeding below the deep chlorophyll maximum. Given the limitations of the selective inhibition method using chloramphenicol, we estimate that fungi contribute no more than 18-22% to total N 2 O production from the oxycline to the oxic-anoxic interface. It remains challenging for omics-based approaches to provide molecular evidence for fungal denitrification, partially because current databases and recent studies have primarily focused on terrestrial fungi. Overcoming the shortfalls of existing methods necessitates approaches that can target fungal diversity and function, such as the use of RNA-seq combined with eukaryotic messenger RNA enrichment or the combination of fluorescence-activated cell sorting and (meta)genome sequencing. These new methods can increase the likelihood of capturing genetic evidence for fungal activities and functional diversity.
Supplementary Materials: The following are available online at https://www.mdpi.com/2309-6 08X/7/3/218/s1, Table S1. Volume of seawater filtered for particulate matter collection. Table S2. Primers designed for PCR. The parts in bold font are standard Illumina adapters. The five ITS3tag primers were combined in equimolar as the forward primer. Table S3. Fungal community composition at the phylum level assessed by the internal transcribed spacer region 2 (ITS2) sequencing. Table S4. (Separate spreadsheet) Taxonomy and read count of each amplicon sequence variant (ASV) in each sample. Table S5. N 2 O production rates integrated from the oxycline to the oxic-anoxic interface depths and the percentage of potential fungal contribution. Table S6. Relative abundance of putative denitrifiers as part of the total fungal community assessed by the internal transcribed region 2 (ITS2). Non-zero values are highlighted in bold. Table S7. (Separate spreadsheet) Scientific names and kingdoms of the top blastp result against NCBI nr database for genes identified as putative P450nor by a hidden Markov model search. All genes were from previously published metagenome assemblies sampled at a station near Station 2 in this study (see Figure S1). Table S8. Comparison of actual number of reads and the estimated minimum number of reads required to recover one copy of P450nor in each sample from the Fuchsman et al. (2017) study. The estimated minimum number of reads required were calculated assuming an average fungal genome size of 40,973,539 bp and the length of the P450nor gene as 1056 bp. Figure S1: (a) Sampling location from March 2012 ("Fuchsman17") where the metagenome samples investigated in this study were collected, and sampling locations in this study ("Station 1" and "Station 2"). Color contour shows oxygen concentrations at 100 m depth from World Ocean Atlas 2013 (March average from 1955-2012). (b) Comparison of potential density (σ θ ), oxygen concentration (O 2 ), and fluorescence measured by Seabird profiling sensors from this study (2018) and from the "Fuchsman17" study (2012). Figure S2. Incubation scheme using a selective inhibition method combined with 15 N tracer technique to determine the fungal contribution to N 2 O production. In each incubation, the final concentration of 15 NO 3 − was 3 µM and the final concentration of 15 NH 4 + was 0.5 µM. Chloramphenicol was added to a final concentration of 87.7 mg L −1 . Figure S3. Calibration curve used to calculate the amount of N 2 O present in each sample based on the total areas under m/z of 44, 45, and 46 ("Area All"). Figure S4. The amount (nmol) and bulk δ 15 N (‰) of N 2 O measured from incubation samples in which N 2 O production rates were non-zero (black empty circles) and of N 2 O concentration standards (blue filled diamonds). Error bars represent standard deviations (n > 5). Figure S5. (a) The fraction of metagenome reads with a DIAMOND hit. (b) The percentage of metagenome reads that were classified as archaea, bacteria, viruses, and algae. The samples were collected during the same month in 2012 at a station close to the station 2 in this study ( Figure S1). Read classification was performed using DIAMOND against the NCBI nr database with an e-value threshold of 1 × 10 −5 .

Acknowledgments:
We are indebted to members of Bess Ward's and Karen Casciotti's Labs, the crew of R/V Sally Ride, and Frank Kinnaman for general assistance. We thank Alyson Santoro's lab and the Analytical Lab at the Marine Science Institute in UCSB for advice on mass spectrometry method development.

Conflicts of Interest:
The authors declare no conflict of interest.

Extrapolation of N 2 O Production by a Putative Fungal Denitrifier
ASV8 was 96.95% similar to the ITS2 region of the N 2 O-producing Penicillium melinii [20], which was isolated from seawater [48]. The relative abundance of ASV8 ranges from 0.1 to 8.5% at hypoxic and anoxic depths (Table S4), and it belongs to the same family (Aspergillaceae) as the N 2 O-producing fungal strain isolated from the Arabian Sea OMZ [21]. We estimated the maximum potential contribution of ASV8 to N 2 O production by scaling its mass-dependent N 2 O production rate (assumed to be identical to Penicillium melinii [20]) by its estimated mass-based abundance, assuming the maximum observed (ITS-based) fractional abundance of 8.5% and fungal carbon content from a similar location. The following equation was used: where P pm is the measured N 2 O production rate for Penicillium melinii [20], m f is the estimated fungal carbon mass (0.12 µg L −1 ) as measured for the eastern tropical South Pacific oxygen minimum zone [10], f ASV8 is the largest fractional contribution of ASV8 to our data set (8.5%), and R N 2 O is the estimated maximum rate of N 2 O production by ASV8. To make R N 2 O an upper bound, we used the largest value of fungal carbon content and an N 2 O production rate by P. melinii (2 mg N 2 O-N g −1 week −1 ) reported previously [20]. By this approach, the R N 2 O extrapolated for this putative fungal denitrifier is 2.1 × 10 −4 nmol L −1 d −1 . This is accounts for at most 0.3% of any fungal N 2 O production rate measured in this study.