Diversity of an uncommon elastic hypersaline microbial mat along a small-scale transect

Sampling site and sample collection

The Archean Domes (26°49′41.7″N, 102°01′28.7″W) is a seasonal, water-fluctuating oval pond (50 m × 25 m (Madrigal-Trejo, 2022)). During the dry season (i.e., most of the year at CCB; Montiel-González et al., 2018), the microbial mats of the Archaean Domes remain humid under a salty crust, which is then dissolved during the rainy season, when the pond fills with water up to ≈20 cm, allowing for the development of the dome-like structures.

The microbial mats of the Archaean Domes were sampled during the rainy season of 2018 (October), with ten evenly spaced samples (S1–S10) along a 1.5-m transect, from the border toward the center, in the ‘Pozas Azules’ ranch of Pronatura Noroeste within CCB (Fig. 1). Salinity and pH were constant along the 1.5-m transect (salinity 5.25%, pH 9.8), and the average nutrient content showed only a slight stoichiometric imbalance (C:N:P 122:42:1), as Medina-Chávez et al. (2019) reported.

Microbial mats were collected and transferred to sterile conical tubes (50 mL), stored at 4 °C, and subsequently frozen in liquid nitrogen until processing in the laboratory. Sampling was under the collection permit SGPA/DGVS/03188/20 issued by Subsecretaría de Gestión para la protección Ambiental, Dirección General de Vida Silvestre.

Total DNA extraction and amplicons sequencing

Prior to DNA extraction, approximately 0.5 g of each mat sample was rinsed with sterile water, removing associated sediments. Total DNA extraction was performed following the protocol reported in De Anda et al. (2018).

DNA was sequenced at the Laboratorio de Servicios Genómicos, LANGEBIO (http://langebio.cinvestav.mx/labsergen/) for the V3 region of the 16S rDNA gene (357F universal primer 5′-CTCCTACGGGAGGCAGCAG-3′/519R universal primer 5′-GWATTACCGCGGCKGCTG-3′; (Lane, 1991)), and the ITS DNA region (ITS1F 5′-TCCGTAGGTGAACCTGCGG- 3′/ITS4R 5′-TCCTCCGCTTATTGATATGC-3′ for ITS; (White et al., 1990)). For each sample, 3 µg of genomic DNA (DO 260/280 1.8) was used; 16S amplicons were sequenced within Illumina Next-Seq 500 platform (2 × 150 PE, using forward reads; reverse reads were discarded from the analysis due to low quality), and the ITS amplicons were sequenced with a Mi-Seq platform (2 × 300 PE).

Bioinformatics analyses

The 16S rDNA amplicon reads were quality filtered and dereplicated using the QIIME 2 bioinformatics platform version 2019.7 (Bolyen et al., 2018). Divisive Amplicon Denoising Algorithm 2 (dada2) denoise-single was used; sequences were demultiplexed and truncated at the 20^th bp from the left, and the 120^th bp from the right. ASVs were then taxonomically classified using the silva_132_99_v3v4_q2_2019-7 database (https://www.arb-silva.de/). For all diversity analyses the final ASV table was rarefied to 760,000 reads per sample.

The ITS amplicon reads were processed using the ITS-specific workflow of the dada2 v1.13.1 package (Callahan et al., 2016) in R (R Core Team, 2020). Briefly, cutadapt was used to remove primer and adapter contamination from the raw reads (Martin, 2011). Prior to determining ASVs, reads with inconsistent bases were excluded, and reads were required to have less than two expected errors based on their quality scores. ASVs classification was inferred from forward and reverse reads using the run-specific error rates.

Due to the wide range of sizes of the ITS amplicons reported for fungal species in CCB (485–782 bp; (Velez et al., 2016)), two different strategies were implemented to define taxonomic diversity. First, the post-denoised ASVs were merged after excluding reads with an overlap of less than 12 bp and reads with mismatches in the overlapping region. Chimeras were removed using the consensus method of “removeBimeraDenovo” implemented in dada2. The taxonomic assignment was performed with the UNITE database (Nilsson et al., 2019), using the RDP naive Bayesian classifier (Wang et al., 2007) available in the dada2 package with a minimum bootstrap value of 80. ASVs with a minimum of 85% identity and with E-values lower than e−45 were retained (Merge ITS data: Table S1). ASVs without taxonomic assignment (mostly non-fungi organisms) were classified using blastn (non-redundant nucleotide collection (nr/nt) database) (Altschul et al., 1990), where only those ASVs within identities greater than 90% were enlisted (Table S2). The remaining ASVs were considered non-classified.

As a great proportion of the reads failed to merge (~56%), a parallel classification was implemented. The RDP naive Bayesian classifier was used on the two sets of reads (forward and reverse) to check for congruence in the classification between the ASVs. ASVs that had an assignment bootstrap value lower than 50, or showed inconsistent classifications among forward and reverse reads, were aligned to the NCBI database (non-redundant nucleotide collection (nr/nt) database), to corroborate their taxonomic assignment. A more robust classification of ASVs was reached when using forward ASVs, thus all subsequent analyses were performed using the Forward-Only ASVs; even though some ASVs reached a genus level assignment, the taxonomic classification of forward reads was limited to the phylum and kingdom level (ITS table on Tables S3 and S4); the Eukaryotic classification followed the UNITE database. The ITS community composition table (Forward-Only data) was rarefied to 164,820 reads using the vegan package (Oksanen et al., 2013) in R (R Core Team, 2020).

In sum, the filtered 16S rDNA amplicon sequences included 760,326–1,973,982 high quality reads per sample, with a total of 6,063 ASVs (Table S5). The filtered ITS amplicon sequences included 164,822–339,394 high quality reads per sample, with a total of 923 ASVs (Table S4). ITS amplicons were obtained for only six samples (S1, S2, S3, S4, S6 and S7) due to low DNA amount that remained in four samples (S5, S8, S9 and S10), after 16S rDNA amplicon sequencing.

Alpha and beta diversity estimates

To evaluate all samples from 16S (S1–S10) and six samples (S1–S4, S6, S7) from the Forward-Only ITS data, we generated rarefied diversity plots using the vegan package (Oksanen et al., 2013) in R (R Core Team, 2020).

In order to visualize and count the number of overlapping ASVs within the ten samples from the 16S rDNA gene and the six samples from the ITS region, we constructed an UpSet plot using the UpSetR package (Conway, Lex & Gehlenborg, 2017) in R (R Core Team, 2020).

Alpha diversity (Shannon index) and beta diversity (Bray-Curtis, Jaccard distances) were estimated using the 16S rDNA ASV community composition table (Table S5), as implemented in QIIME 2. For the Forward-Only ITS ASV community composition table (Table S4), alpha and beta diversity were estimated with the phyloseq package (McMurdie & Holmes, 2013) in R (R Core Team, 2020). The beta diversity distance matrix was used to construct an UPGMA dendrogram of samples (Carteron et al., 2012). Mantel tests (999 permutations) (Legendre & Legendre, 2012) were conducted with the vegan package (Oksanen et al., 2013) to unveil the relationships between beta diversity and geographic distance, as well as between beta diversity and the metabolomic profile distance matrix (see below).

Phylogenetic diversity

We estimated phylogenetic diversity within samples using Faith’s Phylogenetic Diversity index (Faith, 1992) and the Phylogenetic Species Clustering index (Helmus et al., 2007) using the picante package (Kembel et al., 2010) in R (R Core Team, 2020). Given a phylogenetic tree, phylogenetic diversity (PD) estimates the total length of branches encompassing a subset of taxa (e.g., sampled communities), where PD increases as the subset of taxa is more phylogenetically diverse (Faith, 1992). Phylogenetic Species Clustering (PSC) is a measure of the degree of clustering of taxa across the phylogeny (Helmus et al., 2007), where a PSC value approaching zero indicates that the taxa analyzed were phylogenetically clustered.

To estimate the two metrics, we used the rarified diversity matrices and the corresponding phylogenetic trees for bacterial and fungal ASVs, separately. We tested the observed patterns of PSC across samples against 1,000 replicates under a null model and estimated standardized effect sizes (SES) for this index (null model Phylogenetic Species Clustering, nmPSC). This allowed us to identify samples that had more or less phylogenetic clustering than expected by chance alone. To construct the null model, we randomized the entries of the rarified diversity matrices while keeping total sample diversity constant. In addition, we estimated PD and PSC for the most abundant and rare ASVs within each site. To define rare (low proportionated or low read count) and abundant (high proportionated or high read count) ASVs, we estimated the first and third quartiles (1Q, 3Q) of the distribution of ASV (organized by low to high read count), for each site separately. We used the sites’ 1Q and 3Q as the thresholds to define rare and abundant ASVs; the sites’ 1Qs ranged from 3 to 16 reads per ASV, whereas the 3Qs ranged from 37 to 306 reads per ASV; rare and abundant ASV numbers using this method were equivalent per site. We estimated SES for PD and PSC as defined above, using matrices for rare and abundant ASV.

The phylogenetic tree was generated by QIIME 2 (Bolyen et al., 2018) for 16S data. For the ITS data (Forward-Only data), a phylogeny was obtained by applying an equal phylogenetic QIIME 2 methodology, using Maftt (Katoh & Toh, 2010) and FastTree (Price, Dehal & Arkin, 2009) within Cipres Science Gateway (Miller, Pfeiffer & Schwartz, 2010) with default XSEDE parameters.

Metabolomics

Untargeted metabolomics

The HRMS-MS/MS raw data was individually aligned and filtered with MZmine 2.17 software (http://mzmine.sourceforge.net/) (Pluskal et al., 2010). Peak detection was achieved as follows: m/z values were detected within each spectrum above a baseline, a chromatogram was constructed for each of the m/z values that spanned longer than 0.1 min, and finally, deconvolution algorithms were applied to each chromatogram to recognize the individual chromatographic peaks. The parameters were set as follows for peak detection: noise level (absolute value) at 1 × 10⁶, minimum peak duration 0.5 s, tolerance for m/z variation 0.05, and tolerance for m/z intensity variation 20%. Deisotoping, peak list filtering, and retention time alignment algorithm packages were employed to refine peak detection. The join align algorithm compiled a peak table according to the following parameters: the balance between m/z and retention time was set at 10.0 each, m/z tolerance at 0.05, and retention time tolerance size was defined as 2 min. The spectral data matrix (comprised of m/z, retention time, and peak area for each peak) was imported to Excel (Microsoft) for heatmap analysis using R (R Core Team, 2020).

Molecular networking was performed from the mzML files using the standard Global Natural Products Social (GNPS) molecular networking platform workflow with the spectral clustering algorithm (http://gnps.ucsd.edu; Aron et al., 2020). For spectral networks, parent mass of 0.01 Da and fragment ion tolerance of 0.02 Da were considered. For edges construction, a cosine score over 0.70 was fitted with a minimum of four matched peaks, and two nodes at least in the top 10 cosine scores (K) and maximum of 100 connected the components. After that, the network spectra were searched against GNPS spectral libraries, and graphic visualization of molecular networking was performed in Cytoscape 3.8.0 (Shannon et al., 2003). Chemical structural information within the molecular network was obtained using the GNPS MolNetEnhancer workflow (Ernst et al., 2019), which incorporated in silico structure annotations from the GNPS library search (metabolomic_support_info.zip in Supplemental Files). The chemical profiles were also manually dereplicated using UV-absorption maxima and HRMS-MS/MS data against the Dictionary of Natural Products v 29.1 and Dictionary of Marine Natural Products 2019 (Taylor and Francis Group, Abingdon, United Kingdom), MarinLite (University of Canterbury, Christchurch, New Zealand) and SciFinder (CAS) databases as described by El-Elimat et al. (2013), for fungal and bacterial small molecules. For candidate search, exact mass accuracy was set to 5 ppm (Sumner et al., 2007). Heatmap visualization of metabolomic data table was generated in MetaboAnalyst (Pang et al., 2020).

Microbial taxonomic composition

Overall, 40 prokaryote phyla were detected from the 16S rDNA gene libraries (Fig. S1A), and a total of 6,063 different ASVs were found, with a total proportion of 11.8% of unclassified reads (Table S3). Bacteroidetes was the most abundant (high read count) phylum, representing around 23% of all sequences in the samples, followed by Proteobacteria with ≈19% and Cyanobacteria with nearly 18%. Spirochaetes and Chloroflexi represented ≈8% and ≈4%, respectively. Nearly 2% were Firmicutes, Patesibacteria and Halanaerobiaeota (Fig. S1A).

Genera composition is shown in Fig. 2. A total of 219 genera of Bacteria were detected in the 10 sampled sites (S1–S10) at a 1.5-m scale. From those, only nine genera had high read counts (representing ≈80% of the total proportion, Fig. 2A bottom): Coleofasciculus (formerly Microcoleus), Spirochaeta, Catalinimonas, Desulfovermiculus, Halanaerobium, Tangfeifania, Sungkyunkwania, Sediminispirochaeta and Imperialibacter. The remaining 210 genera can be classified as rare biosphere (<1% relative abundance each), with ≈20% of the total proportion (Fig. 2A, top), having a heterogeneous distribution between sites (percentage of relative abundance and percentage of total proportion are used as synonyms).

As for 16S rDNA ASVs, 18 of the 6,063 were abundant (>1% relative abundance each), representing 49% of the relative abundance of reads (Table S5). This 18 ASVs were distributed within the previously mentioned highly abundant (high read count) genera and phyla (except Firmicutes and Patesibacteria).

Conversely, 6,045 ASVs were rare taxa (<1% relative abundance each) which represents 51% of total reads: 4,124 of them were non-classified, and the rest (1,921) were distributed within the whole 40 bacterial phyla and almost all (217) genera (Table S5). The two exceptions of genera with no rare ASVs were Tangfeifania with two highly proportionated ASVs, and Coleofasciculus with only one high-read count ASV. Moreover, 28 ASVs of Archaea were obtained using the 16S primers, comprising two phyla, Woesearchaeota and Euryarchaeota (Fig. S1A). We detected four genera within Halobacteria class: Halovivax, Halorubrum, Halorussus and Haloasta.

For the ITS region libraries, we successfully obtained reads for six of the sampled sites (S1 to S4, S6 and S7). For the Forward-Only ITS data, we detected 923 different ASVs (Table S4). Many of the classified reads belong to Alveolata protists (77%, Fig. S1B), followed by Fungi (11%) and Viridiplantae (8%); Metazoa, Protista, Amoebozoa, Rhodoplantae, Chromista and Stramenophila comprised the remaining 4%. Unassigned taxa were 50% of all reads (Table S3). The majority of ASVs within Forward ITS data could not be accurately classified at lower taxonomic categories, particularly for three common ASVs (ASV0 (unclassified), ASV1 (Ciliophora) and ASV2 (Ciliophora)), which comprise nearly 70% of the proportion of ITS reads for all sites.

Using Merge ITS ASVs data (Table S1) for a better taxonomic assignment, Fungal ASVs comprising 18 genera were detected, with phyla Ascomycota, Basidiomycota and Rozellomycota present in all samples (Fig. 2B). Non-fungal Eukarya Merge ITS sequences classification (hits <90%) are shown in Table S2. Fabrea salina (Alveolata), and Vannella simplex (Amebozoa) were two species that constituted nearly 0.5% of the total proportion of Merge ITS ASVs data, among other classified but low- read count eukaryotic ASVs.

Alpha and beta diversity

Rarefaction curves (Figs. S2A and S2B) showed similar results for both the Forward-Only ITS region libraries (S1, S2, S3, S4, S6, S7), and the 16S rDNA gene libraries, as all sites reached an asymptote, suggesting an adequate sampling effort. Only 183 ASVs (3% of the total ASV read count) were shared between all sites for 16S rDNA gene libraries (Fig. 3A), forming a very small core, but comprised 72.5% of all reads (composed of twenty-one phyla), while most of the ASVs were unique to a specific site or shared between two or three sites. For ITS libraries, the 56 shared ASVs (6%) comprised 95% of all reads (Fig. 3B).

The 16S rDNA ASVs alpha diversity (Table 1) ranged from 2,399 observed ASVs in S2 (Shannon = 7.4, highest Shannon index value), to 544 observed ASVs in S9 (Shannon = 6.1). In the case of the ITS region libraries (Forward-Only ITS data), the highest ASV number for Eukarya was S6 with 383 ASVs and a Shannon diversity of 2.5; S3 had the lowest observed ASVs with 199, but this site had the highest Shannon value of three (Table 1).

Pairwise dissimilarities between samples for both the 16S rDNA gene and ITS region data ranged from 0.2 to 0.7 (16S) and 0.2 to 0.6 (ITS) for the Bray-Curtis distance values and 0.5 to 0.8 (16S) and 0.3 to 0.7 (ITS) according to the Jaccard-1 values (Figs. S3A and S3B). For both distance values, S10 is different from the rest of the 16S rDNA gene libraries, and S1, S2, S4, S6, S7 were associated in the dendogram clusters. In the ITS libraries, there are two main clusters (S1, S2 and S4; S3, S6 and S7). A significant geographic pattern in the Mantel test (p = 0.0402, r = 0.3224) was found for Jaccard-1 of the 16S library, however we did not find significant values based on geography for ITS libraries for either Jaccard-1 or Bray-Curtis (Table S6).

Phylogenetic diversity

For both libraries, Faith Phylogenetic Diversity indicated that rare (low read count) ASVs on each site, were more diverse than abundant (high read count) ASVs (Table 1).

In respect to the null model distribution of Phylogenetic Species Clustering (nmdPSC), the 16S and ITS communities analyzed per site had in some sites a phylogenetic clustering distribution, and in others sites a null model distribution compared with chance alone (Figs. 4A and 4B, dots between lines or below −2).

When analyzing rare and abundant subcommunities using an equivalent ASV number per site, nmdPSC indicated that abundant ASVs were phylogenetically clustered within all sites and for both libraries, while rare ASVs were subcommunities with a null model distribution in some sites, or phylogenetically dispersed in other sites, for bacterial and fungal libraries (Figs. 4A and 4B, Table S7).

Metabolomic profile

Based on the LC-MS results of S1–S10, a wide range of metabolites were visualized on a heatmap analysis chart (Figs. S4A and S4B) and a GNPS molecular networking analysis (Fig. 5 and Table S8).

Untargeted metabolome analysis of microbial communities in soil or environmental samples is difficult, due to the complexity and heterogeneity of this matrix and the scarcity of the secondary metabolites extracted. In the heatmap, the 1,341 molecular features at unique retention times (Figs. S4A and S4B) shared several common regions across all samples (molecular features), allowing us to group them into three main clusters (cluster 1, samples 3–5 and 10; cluster 2, samples 1, 2 and 6; and cluster 3, samples 7–9, Fig. S4B). We found highly heterogeneous metabolic profiles in all samples, and no significant associations were observed between metabolite features and ITS or 16SrDNA community composition matrices (Table S6).

Moreover, feature-based GNPS molecular networking analysis grouped the metabolite features in 1,534 nodes distributed into 67 chemical families (sub-networks) of ≥3 members, 72 families of two metabolite features, and 765 singletons (Fig. 5). Chemical ontology analyses revealed the presence of superclasses of lipids and lipid-like compounds, organic acids, and organonitrogen and organoheterocyclic compounds, all commonly found in soil samples (Fig. 5). Detected compounds in the GNPS and manual metabolome analysis of the main ions observed in LC-MS/MS data focused on microbial metabolites, allowing tentative identification of several secondary metabolites that passed our quality criteria (≤5 ppm mass error; Table S8). Finally, some plant-derived secondary metabolites and small molecules reflecting human activity were also observed (Table S8).

Microbial composition at a small scale

Our study describes microbial diversity on a 1.5-m scale in a rare elastic microbial mat located at an extreme site. Altogether, we found 40 prokaryotic phyla, including 219 genera of Bacteria and 4 of Archaea, as well as 18 Fungal genera (Fig. 2 and Fig. S1).

Alpha and beta diversity

The Shannon Index accounts for both species evenness and abundance; hence, it might not reflect the total richness if most of the found taxa are part of the rare biosphere. For example, Site three for the Forward-Only ITS data had the highest Shannon value, despite having the lowest ASV number (Table 1), due to its higher evenness distribution of taxa. Meanwhile, the rarefaction and observed diversity reflects the richness of the site (Figs. S2A and S2B; Table 1).

Beta diversity analysis showed values that represent a large ASV turnover between sites despite the small scale (Fig. S3), illustrated by the heterogenous composition observed within bacterial genera (Fig. 2A top). Accordingly, the different domains analyzed (i.e., the bacterial 16S rDNA data or the Eukarya Forward only ITS data) had a distinctive association pattern. This turnover pattern should be taken cautiously, as ASVs turnover might a methodological artifact: beta diversity analysis uses data on relative abundances, based on compositional data; this means that the presence of abundant ASVs might “push” the rare ASVs below the detection threshold. Therefore, this pattern should be corroborated in a future time-series using metagenomic data, and additionally, a transcriptome approach would explore if the rare ASVs are an active part of this community. A statistically marginal but significant geographic association was found only for the bacterial 16S rDNA data (p = 0.0402, r = 0.3224); this suggests a geographic pattern that also should be explored with a more extensive sampling.

Bacterial community comparison with other CCB layered microbial communities

We compared our data with four previously published metagenomes of layered microbial communities from the CCB: one microbialite from a large spring-fed river Rio Mesquites (Nitti et al., 2012), one stromatolite from a blue pool in Pozas Azules (Breitbart et al., 2009) and two microbial mats, one from a small and shallow pond called Lagunita in Churince (De Anda et al., 2018), and the other from a green pool in Pozas Azules (Bonilla-Rosso et al., 2012).

At the phylum level, we found seven taxa present in all these five microbial communities: Actinobacteria, Bacteroidetes, Chloroflexi, Cyanobacteria, Firmicutes, Planctomycetes and Proteobacteria. All of them shared among all sites on Archaean Domes mats and had a high read count (Fig. 3A).

Additionally, Bonilla-Rosso et al. (2012) reported 28 orders that were shared within the microbial mat from a green pool (Bonilla-Rosso et al., 2012), the stromatolite (Breitbart et al., 2009) and a coastal hypersaline non-thermophilic microbial mat from Guerrero Negro (Kunin et al., 2008). We found 18 of those 28 orders in the Archaean Domes community (Table S9).

We also explored lower taxonomic levels but did not find a consistent pattern.

The common pattern at higher taxonomic ranks reflects, as Bonilla-Rosso et al. (2012) mentioned, that layered microbial communities assemble in biogeochemical gradients related to defined ecological niches (e.g., photosynthesis, sulfate reduction, heterotrophy) according to their functional traits, rather than their species (Burke et al., 2011, Escalas et al., 2019).

Phylogenetic diversity

The null model distribution of Phylogenetic Species Clustering (nmdPSC, Fig. 4) of whole bacterial and eukaryotic communities showed negative or null values, indicating a clustered pattern compared with the null model. Phylogenetic clustering is predicted to be evident in environments with poor nutrient availability (Mondav et al., 2017), such as the microbial mat conditions. Since little is known about patterns and processes of rare and abundant bacterial and eukaryotic taxa in microbial mats, we also analyzed these two microbial subcommunities.

Metabolomic profile

The dynamic and complex metabolomic profiles (Figs. S4A and S4B, Fig. 5 and Table S8) of the studied microbial mat was consistent with the complex microbial community observed in the metagenomic analysis (Medina-Chávez et al., 2019).

Even though the expression of regulatory genes responsible for secondary metabolites production depends on several factors (e.g., nitrogen or phosphate starvation, stressors like heat, pH and damage to the cell wall, among others) (Baral, Akhgari & Metsä-Ketelä, 2018), along with the lack of soil metabolomics studies and the solubility of the compounds for LCMS analysis, we were still able to tentatively annotate several classes of compounds with interesting biological properties (Fig. 5 and Table S8). Many of these compounds were found to overlap between most of the samples: the antibiotics 8,9-dihydroindanomycin (Li, Roege & Kelly, 2009), neoindanomycin (Rommel, Li & Kelly, 2011) described from Streptomyces species; the cytotoxic daryamide C (Asolkar et al., 2006), gibbestatin B (Kinashi & Sakaguchi, 1984) and furaquinocin H (Ishibashi et al., 1991), from Streptomyces strains; the germination inhibitor colletofragarone A1, from fungi Colletotrichum fragariae and C. capsiciae (Inoue et al., 1996); periconiasin G (Zaghouani et al., 2016) from S. nitrosporeus and Periconia sp. F-31 with free radical scavenger and weak anti-HIV activities; and phomamide and tomaymycin I from Phoma lingam (Ferezou et al., 1980) and Nocardia sp. (Tozuka, Takasugi & Takaya, 1983) with no biological activity described. Other metabolites of interest were found in specific samples, for example the antifungal agents griseofulvin and dechlorogriseofulvin, produced by microorganisms of the genera Penicillium, Aspergillus, Xylaria, Coriolus, Palicourea, Memnoniella, Khauskia, Carpenteles, Arthrinium, Nigrospora, Stachybotrys, and Streptomyces (Atta-ur-Rahmanro, 2005; Kimura et al., 1992; Knowles et al., 2019; Park, Cho & Simpson, 2006) were identified only in S7; the antimitotic agent nostodione A from the blue-green alga Nostoc commune (Kobayashi et al., 1994) in samples S2 and S6–S8; and the mycotoxins citreohybridone B (Kosemura, 2003), setosusin (Fujimoto et al., 1996), roridin E (Böhner et al., 1965) and ganoboninone A (Ma et al., 2015) in samples S8–S10. Interestingly, several observed ions (Fig. 5) do not match with any data reported and could represent potentially new secondary metabolites. Finally, plant- and animal-derived metabolites and compounds resulting from human activity (contaminants), were also annotated (Table S8).

The association of compounds between fungi and osmotrophs, as well as to other Eukarya, reflect the importance of studying the potential role of microorganisms besides bacteria within microbial mats (Carreira et al., 2020). It is important to mention that metabolomics studies of samples containing diverse communities are scarce; thus, this work significantly expands our understanding of soil metabolites, especially from elastic hypersaline microbial mats.