Even though our understanding of marine fungal diversity is increasing [1, 2], a comprehensive knowledge of their active functional ecology remains limited, especially in the open ocean [1]. Fungal activity has been detected in corals, deep sea and coastal sediments, and associated with phytoplankton blooms, including parasites [1, 2], but the full extent that fungi are functionally active throughout the open ocean water column is yet to be established.

Many terrestrial fungi occupy key roles as saprotrophs by decomposing and recycling biogenic matter, making them intrinsic components of healthy functioning ecosystems [3]. In coastal waters there is evidence that planktonic fungi (mycoplankton) degrade and utilise phytoplankton-derived carbohydrate-rich matter in broadly analogous functional modes [4]. However, the extent that carbohydrate-based fungal saprotrophy occurs in the open ocean remains largely speculated [1].

Glycoside hydrolases (GHs) are a widespread group of carbohydrate-active enzymes (CAZymes) [5] that degrade complex polysaccharides and are categorised into substrate-specific families. In terrestrial fungi, secreted CAZymes are key to the functional potential of saprotrophs and are the primary mode of degradation of high-molecular weight (HMW) polysaccharides (e.g. cellulose). Coastal saprotrophic mycoplankton also employ secreted GHs to degrade phytoplankton-derived HMW carbohydrate-based substrates [4], but the prevalence and identity of the specific GH families of active open ocean mycoplankton are unknown.

Metagenomes from the Tara Oceans project have been used to assess mycoplankton diversity [6, 7], but the associated metatranscriptomes are yet to be fully explored from a fungal perspective. We interrogated the Marine Atlas of Tara Ocean Unigenes (MATOU) metatranscriptomic occurrences database [8] for transcribed fungal GH genes to explore broad-scale depth-dependent structuring in the oceans. The MATOU database consists of all unique eukaryotic genes assembled from the Tara Oceans metatranscriptomes (unigenes), their associated taxonomy and occurrence within samples (full methods described in [8]).

To identify GHs within the MATOU database, reference libraries were created for 61 fungal GH families and clustered using CD-Hit [9]. The MATOU unigenes were searched against these libraries using Diamond v.0.9.22 [10], yielding a database where each positive unigene match represents a unique GH (similar genes clustered when similarity < 95% over 90% of the smallest sequence [8]). The database was screened for non-fungal unigenes using the MATOU taxonomy (Fig. 1a, Supplementary Fig. 1). After removal of redundant matches (i.e. where multiple GHs matched to a single unigene), 1,326 unique fungal GH unigenes were found (~0.001% of the entire unigene catalogue) that occurred 44,386 times in all Tara Oceans samples.

Fig. 1: Bioinformatics pipeline and glycoside hydrolase unigene abundance.
figure 1

(a) Pipeline describing the steps involved in identifying fungal CAZymes within the Tara Oceans MATOU database. A fungal GH protein sequence reference database was created from all the 61 characterised GH subfamilies found in fungi. The database was consolidated by clustering sequences at 95 % identity using CD-Hit before Diamond BLAST databases were generated for each subfamily. Unigenes were searched against each of these 61 databases using the following thresholds: e value > 1e−30, score > 1, subject Cov > 75%, keeping only the best alignments. Positive matches were then screened using the MATOU taxonomy to discriminate between fungal and non-fungal unigenes. Occurrences of each unigene within the Tara Oceans transcriptomes were returned. (b) Fungal GH groups found in the MATOU database ranked by abundance over all Tara Oceans samples. (c) Total numbers and taxonomy (including Ascomycota (green), Basdiomycota (orange), and Unassigned (yellow)) of unique fungal unigenes from the ten most abundant GH groups.

Full size image

The top ten GH families containing the greatest number of unique genes were determined by ranking the sum of all unigene occurrences from all samples for each family (Fig. 1b). Overall, the greatest number of unique GHs were involved in cellulose/hemicellulose degradation (GH7). Other substrates of the most abundant GH families included β-glucans (GH17 and GH72), β-glycans (GH5, GH16, GH3), α-glucans (GH13), chitin (GH18), N-/O-glycans (GH47) and xylan (GH43). Fungal GHs were dominated by genes originating from the Ascomycota, except for the abundant GH7s, many of which were unclassified fungal genes (Fig. 1c).

Transcribed fungal GH genes were recovered from 66 stations (Fig. 2). The number of unique GH7s was greatest in the surface and deep chlorophyll maximum (DCM), especially in high productivity areas such as the Mediterranean Sea and the Indian Ocean (Fig. 2b). At stations where concomitant samples were available from the surface, DCM and mesopelagic, a distinct drop in the number of unique GH7s was seen in the mesopelagic. The number of unique unigenes in other abundant GH groups was more heterogeneous between sites (Fig. 2c).

Fig. 2: Depth distribution of fungal glycoside hydrolases in the global oceans.
figure 2

(a) Global map indicating Tara Oceans stations searched for fungal GH unigenes in the surface, deep chlorophyll maximum (DCM) and mesopelagic. (b) Mean unique unigenes/station for each of the major oceanic regions sampled. (c) Depth-dependent partitioning of fungal GH unigenes in the surface (SUR), DCM and mesopelagic (MES).

Full size image

The high prevalence of unique GH7 transcribed genes in surface waters is likely a response to increased ‘fresh’ phytoplankton-derived matter in the photic zone. Cellulose is a key structural component of many phytoplankton cells [11]. Polysaccharides, such as cellulose, also represent a primary source of particulate organic carbon (POC) [12]. Of the enzymes responsible for cellulose degradation, those within the GH7 family are typically the most active, and are important in biomass degradation by terrestrial fungi [13]. The decline in the number of unique GH7 genes below the DCM in the mesopelagic zone in the six sites sampled suggests a shift in mycoplankton functionality due to depletion of readily available phytoplankton-derived carbon sources, and corresponds with the decrease in overall polysaccharide concentrations between surface waters and the mesopelagic zone [14].

Amongst the other major GH groups, the greatest number of unique genes was in the GH17 family, a group of CAZymes that degrade β-glucans. Cladosporium mycoplankton isolated from coastal waters have been shown to secrete GH17 β-1,3-glucosidase when utilising laminarin [4], which is an algal-derived β-glucan that is a major POC component [15]. Also prevalent were GH18 genes, responsible for degrading chitin, another major polysaccharide and important component of zooplankton and some phytoplankton, suggesting that chitin degradation is a functional role of marine mycoplankton, as with fungi in freshwater lake ecosystems [16].

While the number of transcribed GHs is a strong indicator of active carbohydrate metabolism in mycoplankton communities, the identity of these fungi remains uncertain. The extent of GH7 genes with unresolved fungal taxonomy opens interesting questions about the phylogeny of these taxa. The majority of the other GHs were affiliated to the Ascomycota, in line with phylogenetic studies that show the phylum dominates open ocean mycoplankton diversity [6]. However, there is a lack of early-diverging taxa (e.g. Chytridiomycota) within the MATOU database. Since Chytridiomycota parasitism of phytoplankton takes place in the open ocean [17] their absence highlights outstanding gaps in our understanding of the functional ecology of marine fungi.

The vertical flux of POC in the open ocean is an essential feature of the biological carbon pump (BCP), sustaining the oceans capability to sequester carbon [18]. Biogeochemical models of the BCP do not currently consider fungi (e.g. [19]). Given that marine snow is an apparent hotspot for fungi [20], and that here we show differences in fungal GH expression suggesting depth-dependent resource partitioning in relation to POC-associated substrates, the recently proposed ‘mycoflux’ [2] should be considered within a contemporary view of the BCP.