Bacterial community structure and metabolic potential in microbialite-forming mats from South Australian saline lakes.

Microbialites are sedimentary rocks created in association with benthic microorganisms. While they harbour complex microbial communities, Cyanobacteria perform critical roles in sediment stabilisation and accretion. Microbialites have been described from permanent and ephemeral saline lakes in South Australia; however, the microbial communities that generate and inhabit these biogeological structures have not been studied in detail. To address this knowledge gap, we investigated the composition, diversity and metabolic potential of bacterial communities from different microbialite-forming mats and surrounding sediments in five South Australian saline coastal lakes using 16S rRNA gene sequencing and predictive metagenome analyses. While Proteobacteria and Bacteroidetes were the dominant phyla recovered from the mats and sediments, Cyanobacteria were significantly more abundant in the mat samples. Interestingly, at lower taxonomic levels, the mat communities were vastly different across the five lakes. Comparative analysis of putative mat and sediment metagenomes via PICRUSt2 revealed important metabolic pathways driving the process of carbonate precipitation, including cyanobacterial oxygenic photosynthesis, ureolysis and nitrogen fixation. These pathways were highly conserved across the five examined lakes, although they appeared to be performed by distinct groups of bacterial taxa found in each lake. Stress response, quorum sensing and circadian clock were other important pathways predicted by the in silico metagenome analysis. The enrichment of CRISPR/Cas and phage shock associated genes in these cyanobacteria-rich communities suggests that they may be under selective pressure from viral infection. Together, these results highlight that a very stable ecosystem function is maintained by distinctly different communities in microbialite-forming mats in the five South Australian lakes and reinforce the concept that 'who' is in the community is not as critical as their net metabolic capacity.

INTRODUCTION

Microbialites are organosedimentary deposits created in association with complex benthic assemblages (mats) of microorganisms (Burne & Moore, 1987). Examples include stromatolites, thrombolites, dendrolites, leiolites and other microbial deposits, which differ by way of form and fabric (Riding, 2000). Most microbialites are comprised of carbonate, with the degree of mineralisation depending on biotic factors, such as microbial metabolism and production of extrapolymeric substances (EPS), as well as abiotic factors, such as changes in water and atmospheric chemistry, which influence carbonate saturation (Riding, 2011).

While microbialite‐forming mats are comprised of diverse phyla, Cyanobacteria are thought to play a major role in the lithification process. They achieve this through oxygenic photosynthesis, which raises the pH of the microenvironment, encouraging carbonate precipitation, as well as through the production of copious quantities of EPS, which trap sediments and stabilise microbial assemblages (Dupraz et al., 2009, 2011). Other organisms may play a secondary role in microbialite accretion by promoting carbonate precipitation through sulphate reduction, denitrification, ammonification, methane oxidation, or ureolysis (Zhu & Dittrich, 2016).

‘Modern’ microbialites are found in diverse and often extreme environments ranging from geothermal springs (Coman et al., 2015) to fresh (Breitbart et al., 2009; Iniesto et al., 2021), marine (Casaburi et al., 2016; Louyakis et al., 2017) and hypersaline (Allen et al., 2009; DiLoreto et al., 2019; Goh et al., 2009; Ruvindy et al., 2016) aquatic habitats and caves. These environments are typically rich in dissolved minerals, such as carbonate, sodium and magnesium salts and have an alkaline pH, which favours precipitation (Paul et al., 2016; Zhu & Dittrich, 2016). Some of the largest and most distinct living microbialites are the stromatolites of Shark Bay; a hypersaline lagoon located on the coast of Western Australia. These approximately 7000–8000‐year‐old microbialites (Logan et al., 1970) are considered extant analogues of the earliest lifeforms on Earth (White, 2020) and are protected by the United Nations Educational, Scientific and Cultural Organization (UNESCO).

A wealth of ‘omics data has been generated for the stromatolites of Shark Bay (Babilonia et al., 2018; Ruvindy et al., 2016) and other microbialites from Western Australia, including those found in Lake Clifton (Gleeson et al., 2016; Warden et al., 2016) and the hypersaline lakes of Rottnest Island (Mendes Monteiro et al., 2020). These studies have shown that the microbial communities associated with microbialites in these regions are generally dominated by Proteobacteria, Cyanobacteria and/or Bacteroidetes, depending on mat morphology and degree of lithification. Further, a recent metagenomic study revealed that the pustular stromatolite‐forming mats in Shark Bay were enriched in photosynthesis pathways, whereas colloform and smooth stromatolite‐forming mats had an abundance of genes associated with heterotrophic metabolisms, such as sulphate reduction (Babilonia et al., 2018). Microbialites‐forming mats have also been identified in several saline lakes on the coast of South Australia, including Lake Sleaford Mere (Eyre Peninsula), Coorong lagoon, Deep Lake and Lake Hamilton (Rosen et al., 1988; Walter et al., 1973; Warren, 1982). Despite the cultural and ecological significance of these regions, which harbour diverse and endangered wildlife (e.g. Mosley et al., 2019), their microbial communities have never been explored. The aim of this study was to identify the dominant bacterial taxa present in these microbialites and explore their potential metabolic activities driving lithification.

MATERIALS AND METHODS

Site description, sampling and environmental variable measurements

A total of 18 microbial mat and sediment samples (Table S1) were collected from five lakes in South Australia between the 24th and 30th of October, 2013, including Sleaford Mere (SM), Coorong South (CS), Coorong North (CN), Deep Lake (DL) and Lake Hamilton (LH). SM and LH are located on the tip and western side of Eyre Peninsula, respectively, DL is located on the tip of Yorke Peninsula, and CS and CN are situated in Coorong National Park, approximately 5 and 16 km south of the township of Salt Creek, respectively (Figure 1). First descriptions of microbial structures and physical settings of the lakes have been reported for CS and CN (Rosen et al., 1988; Walter et al., 1973). Details on SM and DL can be found in an account by (Warren, 1982).

FIGURE 1

Geographic location of the sampling sites with conductivity gradient overlay (top), and site views and morphologies of sampled microbialite‐forming mats (bottom). The map image was generated in R statistical software using the package leaflet. LH = Lake Hamilton, SM = Sleaford Mere, DL = Deep Lake, CN = Coorong North, CS = Coorong South

Submerged microbial mats as well as structures on the water's edge and several centimetres away from it were sampled. No sample was taken from a depth greater than 30 cm to ensure solar irradiation was similar across all samples. All mat samples showed obvious signs of microbial colonisation such as coloured laminations or green subsurface layers (Table S1; Figure 1). Lithified mats and soft mats were further classified based on their mesostructured (internal fabric visible to the naked eye; Dupraz et al., 2011) in which lithified mats showed various degrees of lithification (rock‐like layering), whereas soft mats lacked lithification and were soft and slimy. Details of the mat morphologies are provided in Table S1. Sediment controls consisted of bare submerged sand, devoid of any visible microbial growth.

Water samples were also collected in the vicinity of the microbial structures and conductivity was measured with a KCl‐standardised Horiba 9382‐10d probe. Following the classification of wetlands and deep‐water habitats proposed by Cowardin and colleagues (Cowardin et al., 1985), we categorised the five lakes based on their specific conductivity as mesosaline (SM, 22 mS/cm), polysaline (CS, 37 mS/cm) and hypersaline (CN, DL and LH, 69, 120 and 141 mS/cm, respectively). Unfiltered water pH was measured using a freshly calibrated Dynamica pH MasterBIO (Scientifix, Victoria, Australia) probe. Major cations (Na+, K+, Ca2+ and Mg2+) were measured by ICP‐OES (Agilent 720, Agilent, California, USA) using a ‘OneNeb’ nebuliser using standard operating parameters. Instrument drift was accounted for by frequent recalibration using NIST‐derived, high purity, multi‐element standards (Merck Millipore, Darmstadt, Germany). Trace cations (Al3+, Mn2+, Fe+, Ni2+, Cu2+, Zn2+ and Pb2+) were analysed by ICP‐MS (Agilent 7700x), which was operated in helium mode (3.5 mL min−1) to remove potential mass interferences. In addition to frequent recalibration, instrument drift and matrix interferences were monitored by the in‐line addition of an internal standard containing 45Sc, 72Ge and 89Y. Internal standard recovery was typically 90 to 110%, and samples outside the ±20% range were diluted and reanalysed.

DNA extraction and sequencing

For each sample, approximately 500 mg of material was coarsely ground, and DNA was extracted using the MP Bio FastDNA Spin Kit for Soil (MP Bio, New South Wales, Australia). To increase DNA yield at the final elution step, spin columns were incubated at 55°C and the elution buffer was passed through the column twice. DNA concentration and quality were measured with a NanoDrop TM1000 Spectrophotometer (Thermo Fisher Scientific, Massachusetts, USA) and checked visually under UV irradiation following electrophoresis through 1% agarose gels and staining with ethidium bromide. Barcoded primers, 27F (Lane, 1991) and 519R (Lane et al., 1985), containing 8 nucleotide multiplex identifier tags, were used to amplify the hypervariable V1‐3 region of the 16S rRNA gene. PCR reactions were performed in a final volume of 20 µl, containing approximately 10 ng of DNA, 50 pmol each of forward and reverse primer, 1 mm MgCl2, 1 mm BSA, 1.5 mm dNTPs, 0.3 U BioTaq DNA Polymerase (Bioline, New South Wales, Australia) and 0.1 U Pfu proofreading DNA Polymerase (Promega, Wisconsin, USA) in 1′ Tris‐HCl buffer (Bioline). Thermocycling conditions were as follows: initial denaturation at 95°C for 5 min, 30 cycles of elongation at 95°C for 15 s, annealing at 55°C for 30 s, elongation at 68°C for 60 s and a final elongation step at 72°C for 7 min. Final reactions were treated with 0.05 U of exonuclease I and 0.25 U of shrimp alkaline phosphatase (New England Biolabs, Massachusetts, USA) for 15 min at 37°C to remove residual primers and dNTPs. Enzymes were subsequently deactivated for 20 min at 80°C. Normalisation and pooling of amplicons before sequencing were performed using the SequalPrep Normalisation Plate (Thermo Fisher Scientific). The concentration and quality of PCR amplicons were checked with the Qubit® dsDNA HS assay kit and fluorometer (Thermo Fisher Scientific) and their size was measured using the Agilent 2200 TapeStation (Agilent, USA). The Agencourt AMPure XP bead clean‐up kit (Beckman Coulter, California, USA) was used on the pool to reduce primer dimers. Sequencing was performed at the Ramaciotti Centre for Genomics, at the University of New South Wales, Australia, on the Illumina MiSeq instrument using a MiSeq Reagent Kit v3 (Illumina, San Diego, California, USA) with a 2 × 300 bp run format.

Bioinformatic and statistical analyses

The 16S rRNA gene sequencing data were demultiplexed via the Illumina's BaseSpace cloud computing environment and are available on the NCBI Small Read Archive (SRA) under project accession number PRJNA750862 (ncbi.nlm.nih.gov/bioproject/PRJNA750862/). Paired‐end FASTQ files were pre‐processed, quality‐filtered and analysed using the Quantitative Insights into Microbial Ecology (QIIME) 2 pipeline v2020.2 (Bolyen et al., 2019). Firstly, the DADA2 software package (Callahan et al., 2016) wrapped in QIIME2 was used for quality filtering (median PHRED score ≥25), denoising, joining paired ends and removing chimeric sequences, which resulted in a total of 1,363,396 high‐quality, merged sequences with an average length of 480 bp across 18 samples. The DADA2 denoising algorithm inferred those sequences into 15,517 exact amplicon sequence variants (ASVs), of which 15,268 ASVs with minimum frequency of two were retained for further analysis.

A Naïve Bayes classifier was trained using the QIIME2’s q2‐feature‐classifier plugin (Bokulich et al., 2018) and taxonomy assigned using reference sequences from the Greengenes 13_8 99% OTU full‐length dataset (DeSantis et al., 2006). It has been shown that when a Naïve Bayes classifier is trained exclusively with sequences representing the target region, the accuracy of taxonomic classification of 16S rRNA gene sequences improves (Werner et al., 2012); therefore, we used the 27F/519R primer pair to trim the Greengenes reference sequences to the V1‐V3 region and used this fragment for classifier training. All ASVs that were identified as mitochondria or chloroplast sequences were removed from the feature table. De novo phylogeny was constructed from the remaining 11,114 ASVs (hereafter observed ASVs) by MAFFT alignment v7.475 (Katoh et al., 2002) and FastTree v2.1.10 (Price et al., 2010) as implemented in QIIME2 with default masking options. Rarefaction analysis was performed and visualised in QIIME2 using the alpha‐rarefaction function, which showed that the diversity present in all samples was adequately captured (Figure S1). QIIME2 mapping files (*.qza files) are provided as Files S1 and S2.

To examine the differences in community structure across locations, we performed Principal Coordinate Analysis (PCoA) based on the Bray‐Curtis dissimilarity matrix derived from rarefied, square‐root transformed ASVs. The differences were assessed by conducting a permutational multivariate analysis of variance (PERMANOVA) using the command ‘adonis2’ in the package vegan (Oksanen et al., 2020) in R using the grouping factor of either ‘Lake’ or ‘Conductivity’. Within‐lake heterogeneity of dispersion was assessed by a test of multivariate dispersion (Anderson et al., 2006) using the command ‘permdisp’ in the package vegan (Oksanen et al., 2020). Pairwise testing between locations was conducted using the command ‘pairwise.perm.manova’ from the package RVAideMemoire (Hervé, 2019). To account for multiple testing, p‐values were corrected using the Benjamini‐Hochberg procedure.

To gain insight into functional variations of the bacterial communities and identify putative metabolic pathways important for microbialite formation, taxon‐based metabolic profiling using the PICRUSt2 pipeline (Douglas et al., 2020) implemented in QIIME2 was performed with default settings including EPA‐NG as placement tool to place sequences into reference tree, maximum parsimony as hidden‐state prediction method, and a cut‐off of 2.0 for nearest‐sequenced taxon index (NSTI) values. PICRUSt2 uses ASV/OTU taxonomic classifications and the reference genomes to compute a ‘predictive metagenome’ of the target ecosystem from which community functional capacity can be inferred (Langille et al., 2013). Functional predictions were made using the Kyoto Encyclopedia of Genes and Genomes (KEGG) classification (Kanehisa & Goto, 2000) and the MetaCyc pathway database (Caspi et al., 2016). To examine the differences in metabolic potential across samples, PCoA was performed using either Jaccard or Bray‐Curtis distance matrix generated from rarefied, square‐root transformed predicted KEGG genes. To get stratified outputs showing contributed ASVs and their contribution proportions to each KEGG gene, the ‘metagenome_pipeline.py’ python script in PICRUSt2 (v2.3.0 beta) was used. The stratified table is provided in File S3 and can be viewed in R. The commands used to run the QIIME2 and PICRUSt2 pipelines is provided in File S4.

To statistically quantify the difference in relative abundance of the pathways/functions of interest between microbialites and sediments, we used pairwise Wald tests as implemented in the DESeq2 package (Love et al., 2014). Over‐abundant genes/pathways with Benjamini‐and‐Hochberg corrected p‐values < .05 and log2 fold change >1.5 or < −1.5 were retained for further analysis. Calculation and visualisation of shared ASVs, KEGG genes and MetaCyc pathways were performed using Euler (Larsson, 2018) and ggVennDiagram packages (Gao et al., 2021) in R.

To improve the poor taxonomic resolution for Cyanobacteria obtained with the Greengenes database (gg_13_5), cyanobacterial ASVs were filtered from the QIIME2 feature table and then aligned to reference sequences from the Cydrasil cyanobacteria database (v2.0, https://github.com/FGPLab/cydrasil) using the short‐read alignment algorithm, PaPaRa (Berger & Stamatakis, 2011), placed into the reference tree using the RaxMLB algorithm (Stamatakis, 2014) and visualised with iTOL3 server (Letunic & Bork, 2016). The NCBI database was also used to examine the similarity of those ASVs compared to publicly available cyanobacterial sequences. ASVs that were placed in a branch of the Cyrdrasil reference tree with >70% certainty (Roush & Garcia‐Pichel, 2020) and >94.5% similarity to an identified sequence (Yarza et al., 2014) were considered to represent the same genus.

RESULTS AND DISCUSSION

Bacterial community structure

All five South Australian saline lakes harboured diverse bacterial communities. A total of 48 phyla (8093 observed ASVs) were recovered from the microbial mats, while 44 phyla (4330 observed ASVs) were recovered from the sediments (Table S2). Most of the recovered phyla were common to mat and sediment samples; however, the relative abundance of ASVs assigned to each phylum was distinct (Figure 2, Table S2). The microbial mats were dominated by Proteobacteria (32.2%), Cyanobacteria (26.3%) and Bacteroidetes (21.7%), with Actinobacteria (5.5%), Chloroflexi (3.8%) and Planctomycetes (3.4%) comprising a minor component of the combined mat communities. The sediment samples were also dominated by Proteobacteria (40.5%) and Bacteroidetes (25.7%), but Cyanobacteria only contributed a modest proportion (6.1%) in these communities (Figure 2, Table S2).

FIGURE 2

Bacterial community structure in South Australian saline coastal lakes. Relative abundances of the major taxa contributing to the bacterial communities in mats and sediments, in order of increasing conductivity of lakes. SM = Sleaford Mere, CS = Coorong South, CN = Coorong North, DL = Deep Lake, LH = Lake Hamilton. Sediment samples SA4, SA19, SA7 and SA13 belong to lakes SM, CS, CN and LH, respectively. The four most abundant classes of the phylum Proteobacteria are presented. The low abundance phyla in each group with median relative abundance <2% were grouped as Other

Proteobacteria (Alpha‐, Delta‐ and Gamma‐), Cyanobacteria and Bacteroidetes are a common feature of lithifying and non‐lithifying microbial mat systems from a range of fresh to hypersaline habitats. In the present study, we observed a high and relatively consistent abundance of these core bacterial phyla in four of the five South Australian lakes studied, the exception being the mesosaline lake, SM, which had a relatively low proportion of Bacteroidetes (Figure 2). While the mat samples harboured a similar cohort of phyla, they were very heterogeneous at the lowest taxonomic level with only 4 ASVs shared among the five lake communities (Figure 3; Table S3) and were distinctly clustered according to their location (PERMANOVA; p = .0001, F = 2.0738) and/or conductivity (PERMANOVA; p = .0002, F = 2.2922) on PCoA ordination (Figure 4). These observations are largely in agreement with previous studies of microbialite‐forming mats from hypersaline Shark Bay, Western Australia (Ruvindy et al., 2016), thrombolites from hypersaline Lake Clifton, Western Australia (Warden et al., 2016), thrombolites (Mobberley et al., 2013), stromatolites from Little Darby Island, The Bahamas (Casaburi et al., 2016), microbialites and non‐ lithifying mats from freshwater habitats along the Trans‐Mexican volcanic belt (Iniesto et al., 2021) and stromatolites from the peritidal eastern coast of South Africa (exposed to both fresh and marine water) (Waterworth et al., 2021). Interestingly, sediment samples were also clustered with mat samples according to lake/conductivity.

FIGURE 3

Venn diagram showing the number of shared ASVs among mat samples from five lakes; gradient grey colour filling exhibits the differences in ASV counts among various regions

FIGURE 4

Principal coordinate analysis (PCoA) using Bray‐Curtis distance matrix generated from ASV profiles of mat samples. SM = Sleaford Mere, CS = Coorong South, CN = Coorong North, DL = Deep Lake, LH = Lake Hamilton

The Proteobacteria community in the South Australian microbialite‐forming mats was comprised mostly of Alphaproteobacteria (22%), Gammaproteobacteria (6.8%) and Deltaproteobacteria (3%). In the corresponding sediments, the contributions of these classes were 16.7%, 14% and 8.6%, respectively (Table S4). The most abundant Alphaproteobacteria in the mats were anoxygenic phototrophs, specifically purple non‐sulphur bacteria belonging to the orders Rhodobacterales, Rhodospirillales and Rhizobiales. It has been suggested that these taxa are important for microbialite formation due to their anoxygenic photosynthetic activities that promote carbonate precipitation (Dupraz & Visscher, 2005; Gérard et al., 2018) and their ability to fix nitrogen, even in the presence of nitrogen‐fixing cyanobacterial taxa (Havemann & Foster, 2008). However, as these Alphaproteobacteria were also found in high abundance in the surrounding sediments (Table S4), their involvement in the construction of microbialites in this scenario should be interpreted with caution.

Surprisingly, Deltaproteobacteria were poorly represented in the mat samples compared to the sediments. Members of the Deltaproteobacteria, including dissimilatory sulphate‐reducers, such as Desulfobacterales, are thought to drive the alkalinity engine via the oxidisation of organic carbon and the degradation cyanobacterial EPS, thereby promoting carbonate precipitation (Baumgartner et al., 2006; Visscher et al., 2000). However, the Desulfobacterales accounted for only 0.7% of the combined microbial mat communities. In contrast, Desulfobacterales was the third most abundant order detected in the sediments, accounting for 4.9% of all taxa (Table S5).

The Bacteroidetes community in the South Australian microbialite‐forming mats was mostly comprised of Rhodothermi (12.3%), Cytophagia (5.1%) and Flavobacteriia (2.6%) (Table S4). Rhodothermi possess a high‐salinity niche preference (Zhang et al., 2019), which is consistent with the meso‐to‐hypersaline environment of the five lakes. This order was found in similar abundance in the sediments (11.3%). Flavobacteria and to a lesser extent Cytophagia, can degrade complex organic matter including EPS (Elifantz et al., 2005; Fernández‐Gómez et al., 2013; Zhang et al., 2015), providing metabolic substrates for other heterotrophic organisms as well as a structural backbone for biofilm formation (Nitti et al., 2012). The degradation of EPS also creates localised conditions suitable for carbonate precipitation in the presence of carbonate species (Breitbart et al., 2009; Zhu & Dittrich, 2016). However, a direct role of these taxa in carbonate precipitation and thus microbialite formation could not be established here because Flavobacteriia and Cytophagia were also detected at similar abundance (3.9% and 3.7%, respectively) in the sediments (Table S4).

As mentioned above, Cyanobacteria made up a considerably high proportion of the mat community (26.3%) compared to the sediment community (6.1%). DESe2 analysis at the phylum level also revealed that Cyanobacteria was the only phylum that were significantly enriched in lithified mats compared to sediments (adjusted p = .0002; Table S6). These results were expected as Cyanobacteria play a key role in the accretion and stabilisation of microbialites. The Cyanobacteria community in the mats was comprised mostly of the orders Oscillatoriales, Synechococcales, Chroococcales, Pleurocapsales and Nostocales (Table S5). These orders have previously been identified as important microbialite and mat builders (Águila et al., 2021; Babilonia et al., 2018; Emmanuelle Gérard et al., 2013). At lower taxonomic ranks, the cyanobacterial communities were highly heterogenous across the studied lakes. At the genus level, only Phormidium was common to mat samples from all five lakes. This filamentous genus was very low in abundance, accounting for only 0.15% of all 60 detected cyanobacterial genera (Table S7) and thus was unlikely to play a major role in microbialite construction in the studied ecosystems. A closer look at the most abundant cyanobacterial genera revealed some interesting patterns, in which several groups of genera were found to be lake‐ and/or conductivity‐specific (Figure 5). For example, the mesosaline lake, SM, was enriched in Oculatella, Halomicronema and an unclassified genus belonged to the Chroococcidiopsidaceae family. Caldora and Calothrix were among the genera specifically found in two hypersaline Coorong lakes. At the highest conductivity, DL and LH were enriched in halotolerant cyanobacteria belonging to the genera Euhalothece, Cyanothece, Dactylococcopsis and Halothece. While some genera were abundant in three or more lakes (e.g. Leptolyngbya), their comprising species were lake‐specific.

FIGURE 5

Abundances of the major genera contributing to the cyanobacteria community of mat samples from the South Australian saline lakes, in order of increasing conductivity. SM = Sleaford Mere, CS = Coorong South, CN = Coorong North, DL = Deep Lake, LH = Lake Hamilton. The abundance of ASVs were standardised to the median sequencing depth using the standf function in R. Prefix ‘Un‐’ stands for ‘Unclassified’. The shading represents the enrichment of the genera according to conductivity; blue: mesosaline, yellow: polysaline to hypersaline and orange: hypersaline

These observations imply that microbialite production is achieved by distinct Cyanobacteria communities in the five lakes and suggests a very high degree of functional redundancy among community members. Similar trends have been observed in other ecological settings, in which clear, habitat‐ or host‐associated, taxonomic core communities could not be identified at the OTU level (Burke et al., 2011; Huse et al., 2012), even though in some cases, an extensive degree (70%) of consistency amongst functional traits was observed (Burke et al., 2011). These studies thus support the theory that the assembly of bacterial communities is best explained in term of gene content rather than species content (Burke et al., 2011).

As these results collectively show, no particular taxon is diagnostic of microbialite‐forming mats and given the lack of sediment controls in most previous studies, the association of individual taxa with the construction of microbialites should be interpreted with caution.

Putative metabolic functions of bacterial communities based on recovered taxa

The PICRUSt2 algorithm was employed to predict the metabolic potential of the microbialite‐associated bacterial communities and their corresponding roles in the bioprecipitation of carbonate. Although this is an indirect method for estimating microbial metabolic functions and therefore, has some limitations (Douglas et al., 2020; Langille et al., 2013; Sun et al., 2020), it has been demonstrated to accurately predict the functional complexity of microbial communities across a wide range of ecosystems, including microbialite‐forming mats and other microbial mats from geothermal springs (Coman et al., 2015), freshwater (Iniesto et al., 2021; Ramoneda et al., 2021), marine and hypersaline environments (Casaburi et al., 2016; DiLoreto et al., 2019; Louyakis et al., 2017; Saona et al., 2021).

A total of 11,093 ASVs were used as input for the PICRUSt2 pipeline, of which 10,930 ASVs (accounting for 99.75% of total amplicon reads) with NSTI values below 2.0 were used for prediction (Table S8), indicating that the reference genomes for microbialite bacteria were very well presented in the PICRUSt2 reference tree. Across all lake samples, PICRUSt2 predicted a total of 10,439 KEGG functions and 489 MetaCyc pathways. Surprisingly, the mat and sediment groups shared 10,324 KEGG genes and 483 MetaCyc pathways, which equated to around 99% similarity in the metabolic potential between the two groups (even though they shared only 12% of ASVs; Figure 6a). Consistently, PCoA ordination based on binary Jaccard distance that measures the presence/absence of metabolic components showed no separation between sediment, lithified mat and soft mat samples (Figure 6b). However, PCoA on Bray‐Curtis distance clustered four sediment samples together and delineated them from the microbial mat group (Figure 6c), demonstrating that the abundance of each metabolic component was the determinant factor that differentiated these two groups. Interestingly, the two soft mat samples were distinct from the lithified mat group along the axis PC2, suggesting a transition in metabolic activity occurs as mats become lithified. To obtain a clear distinction between lithified microbial mats and sediments, the soft mat samples were excluded from subsequent analyses.

FIGURE 6

PICRUSt2‐predicted metabolic functions of lithified mats, soft mats and sediments. (a) Euler diagrams represent the number of shared ASVs, shared predicted KEGG‐annotated genes and shared predicted MetaCyc pathways between sediments and the mat group. (b) Principal coordinate analysis (PCoA) using Jaccard distance matrix generated from predicted KEGG gene profiles. (c) PCoA using Bray‐Curtis distance matrix generated from predicted KEGG gene profiles. The circles represent groups with significant pairwise PERMANOVA test: adjusted p = .003 for lithified mat group vs sediment group. (d) Venn diagrams represent the number of shared predicted KEGG‐annotated genes (left) and shared predicted MetaCyc pathways (right) among five PICRUSt2‐predicted microbialite‐associated metagenomes from the five examined lakes; gradient grey colour filling exhibits the differences in KEGG/pathway counts among various regions. SM = Sleaford Mere, CS = Coorong South, CN = Coorong North, DL = Deep Lake, LH = Lake Hamilton. (e) The abundance of selected KEGG‐annotated metabolic genes that were significantly more abundant in lithified mats compared to sediment controls based on DeSEq2 analysis. A full description of gene names is provided in Table S10

Comparison of metabolic components of the lithified mats revealed that 483 (99%) MetaCyc pathways and 10,309 (99%) KEGG functions were shared among the mat communities from the five lakes (Figure 6d). This suggests that key functional guilds are highly stable across a gradient of conductivity (meso to hypersaline), despite community heterogeneity in the lithified mats at the ASV level (Figures 3 and 4).

When the pooled lithified mat data was compared to the pooled sediment data, 26 MetaCyc pathways and 864 KEGG functions were observed as being significantly different between the two groups, as revealed by differential abundance (DESeq2) analysis (Tables S9 and S10). Notably, photosystems I (psa) and II (psb) genes were significantly more abundant (3.6 to 7.7‐fold) in the lithified mats compared to the sediments (Figure 6e; Table S11). This result is consistent with photosynthesis being a major driver of biological‐induced carbonate precipitation during microbialite accretion (Dupraz et al., 2009, 2011; Zhu & Dittrich, 2016). The lithified mats were also highly enriched in genes encoding a range of pigments and antennae proteins, including chlorophyll, allophycocyanin, phycocyanin, phycoerythrin and phycoerythrocyanin, and those involved in photosynthetic electron transport and cytochrome b6/f complex (Figure 6e; Table S11).

Similarly, the near complete family of genes encoding carbon‐concentrating mechanism (CCM) proteins (ccmK/L/M/N/O), bicarbonate transporters (cmpA/B/C) and cyanophycin biosynthesis enzymes (cphA/B) were significantly more abundant in the lithified mats compared to the sediments (Figure 6e; Table S11). The CCM allows cyanobacteria to enrich the amount of CO2 at the site of Rubisco by up to 1000‐fold compared to that in the surrounding medium, thereby elevating pH on the cell surface (Miller & Colman, 1980; Price et al., 2008). This results in an alkaline microenvironment that is conducive to carbonate precipitation (Dupraz et al., 2009). These results suggest that photosynthetic performance was coupled with, and hence strengthened by CCM in the lithified mats.

Consistent with the high abundance of anoxygenic phototrophic Proteobacteria, including purple non‐sulfur bacteria (Alphaproteobacteria), purple sulfur bacteria (Gammaproteobacteria), and to a lesser extent, Acidobacteria and Chloroflexi, in all samples (Figure 2; Tables S2 and S4), many genes involved in anoxygenic photosynthesis (e.g. pufA/B/C/L/M/X, puhA, bchC/E/F/J/O/X/Y/Z) were well‐represented in both lithified mats and sediments. However, the abundance of most of these genes was not statically different between the two groups (data not shown). These observations collectively suggest that oxygenic photosynthesis by cyanobacteria is likely the most critical contributor to carbonatogenesis in microbialites in the South Australian lakes examined.

Cyanobacteria moderate pigment production and transport processes in response to abiotic stress, such as nutrient limitation, osmotic stress and high light (Singh & Montgomery, 2013; Tamary et al., 2012; Yang et al., 2020). Therefore, the relatively high abundance of pigment and antenna protein genes in the lithified mat samples could explain how the resident cyanobacteria cope with daily exposure to high ultraviolet radiation and often extreme osmotic pressure. In line with this, the lithified mat group was enriched in genes involved in two‐component systems of the OmpR family, such as nblS and nblR, which regulate the photosynthetic apparatus in response to high light and nutrient stress (Salinas et al., 2007; van Waasbergen et al., 2002), and manS and manR, which play an important role in Mn2+ homeostasis (an essential component of photosynthetic machinery) when Mn2+ is limiting (Yamaguchi et al., 2002) (Figure 4e; Table S9), as was the case in the five lakes (Table S12). The nblR, manS and manR genes were assigned exclusively to cyanobacterial taxa, the majority of which belonged to the families Cyanothecaceae, Pseudanabaenaceae and Xenococcaceae, whereas nblS was associated with almost all major phyla and was particularly well‐represented in Bacteroidetes and Proteobacteria (the classes Alphaproteobacteria and Deltaproteobacteria). Other two‐component‐system related genes that were highly represented in the lithified mats included those associated with quorum sensing (LuxR family), chemotaxis and twitching motility, which all play an important role in adhesion and biofilm formation (O'Toole & Kolter, 1998). Cyanobacterial circadian clock genes were also enriched in lithified mat samples, including circadian oscillator (kaiA/B/C) and circadian input kinase (cikA) genes, and the Synechococcus adaptive sensor gene (sasA) (Figure 6e; Table S11), which were assigned to cyanobacterial genera with known nitrogen‐fixing representatives such as Cyanothece, Chroococcidiopsis, Gloeocapsa, Leptolyngbya, Microcoleus, Nostoc, Pseudanabaena, Synechococcus and Synechocystis. These genes are thought to allow diazotrophic cyanobacteria to alternate photosynthesis and nitrogen fixation, the metabolisms important for carbonate precipitation (Breitbart et al., 2009).

The high abundance of stress response and sensing genes in the lithified mats suggests a strong capacity of microbialite‐forming bacteria to sense and respond to environmental change. Our results reflect previous studies that identified abundant levels of these genes in microbialite‐forming mat metagenomes recovered from freshwater (Breitbart et al., 2009; White et al., 2015), marine (Casaburi et al., 2016) and hypersaline environments (Charlesworth et al., 2019; Warden et al., 2016). Collectively, these studies suggest that signal transduction pathways are highly conserved at the microbial ecosystem level across different geological settings.

In addition to photosynthesis, other metabolisms commonly associated with carbonate precipitation detected by PICRUSt2 included ureolysis, nitrogen fixation, denitrification, sulphate reduction and methanogenesis. However, only genes involved in ureolysis and nitrogen fixation were significantly more abundant in the lithified mats compared to the sediments. For example, genes encoding the three‐subunit ureases (ureA/B/C), urease‐specific accessory proteins (ureD/E/F/G/H/J) and urea transport system proteins (urtA/B/C/D/E) were enriched in the mats (Figure 4e; Table S11). Most of these genes were assigned to the alphaproteobacterial families Hyphomicrobiaceae, Rhodobacteraceae and Rhodospirillaceace, the cyanobacterial families Cyanothecaceae, Pseudanabaenaceae and Xenococcaceae and the two verrucomicrobia families Puniceicoccaceae and Verrucomicrobiaceae. Ureolytic bacteria possess strong carbonatogenesis capabilities because they cause rapid, widespread increases in the alkalinity of the microenvironment due to the release of ammonia and subsequent production of OH− through ureolysis, which in the presence of calcium, results in high rates of calcium carbonate precipitation (Mitchell et al., 2010; Reeksting et al., 2020). Ureolysis has been identified as the greatest contributor to microbial‐induced carbonate precipitation in cave and hypersaline environments (Zhu & Dittrich, 2016); however, the role of ureolysis in microbialite formation has not been widely studied. To the best of our knowledge, urea metabolism genes have only been previously identified in freshwater microbialite‐forming mats from Pavilion Lake, Canada (White et al., 2015). Therefore, ureolysis can be considered a relatively novel feature of the lithified mats in the South Australian saline lakes.

Genes related to nitrogen fixation (nifD/E/H/K/T/W/X/Z), ferredoxin‐dependent assimilatory nitrate reduction (nirA, narB), nitrate transport (ABC transporter, nrtA, nasF, cynA) and bidirectional hydrogenases (hoxE/F/H/U/Y) were enriched in the lithified mats compared to the sediments (Figure 6e; Table S11). The reduction of nitrogen/nitrate/nitrite can lead to calcite precipitation and thus in some cases can supplement ureolysis to improve calcification efficiency (Zhu et al., 2019). Denitrification associated genes such as nir, nos and nor genes, were found in the mat and sediment predictive metagenomes, but they were either enriched in sediment samples or their abundance was not significantly different between the two groups, suggesting that denitrification does not represent primary sinks for nitrogen within the studied microbialite‐associated bacterial communities.

Finally, lithified mats were highly enriched in CRISPR/Cas system and phage shock genes, particularly cyanobacterial genes (e.g. csc1/2/3, csx3/10; Figure 6e; Table S11), suggesting that microbialite communities in South Australian lakes may undergo selective pressure from viral infection. Recent studies have revealed diverse assemblages of single‐stranded DNA viruses in mats associated with modern stromatolites and thrombolites (Desnues et al., 2008; White et al., 2018) and it has been suggested that viral communities may influence the transition from soft mats to lithified structures by, for example, altering cyanobacterial metabolism through increasing primary photosynthetic production and altering the alkalinity engine and/or EPS towards carbonate precipitation, although the exact mechanisms remain unknown (White et al., 2021). Other defence genes identified in the lithified mats and sediments were those related to antibiotic (e.g. bacitracin, methicillin, tetracycline), heavy‐metal (e.g. arsenate, copper, gold, mercury), multidrug and multiple‐antibiotic resistance. These resistance mechanisms are prevalent in microbialite systems globally (Ruvindy et al., 2016; White et al., 2015). However, given the similar abundance of these genes in mats and sediments, it is likely that heavy‐metal and/or antibiotic resistance is a typical feature of microorganisms inhabiting the high‐salinity South Australian lakes studied. Like many other hypersaline environments, the five studied lakes were enriched in heavy metals, particularly aluminium and iron, of which the concentrations were between 64 and 131 times, and 8 and 215 times higher than those of average seawater, respectively (Table S12; Turekian, 1968). Our results reflect previous studies that have observed metal tolerance in halotolerant and halophilic bacteria (Voica et al., 2016).

CONCLUDING REMARKS

Overall, our results suggest that the formation of microbialites in the five South Australian saline lakes is largely determined through cyanobacterial oxygenic photosynthesis, and to a lesser extent, ureolysis and nitrogen fixation. These key functional guilds remained highly stable in microbialite communities across a gradient of conductivity, namely from mesosaline (SM) to polysaline (CS) and hypersaline (CN, DL and LH) environments, despite the dynamic community composition in these ecosystems. This finding also indicates that the core functions relevant to microbialite formation are not restricted to certain taxonomic groups, but members from different taxa could form shared functional guilds.

The high abundance of stress response and sensing functions including circadian clock genes was a notable feature of the microbialite‐forming mat communities. Given the extreme similarity in metabolic capacities between mat and sediment samples, such genes/functions can be regarded as having a critical role in fine‐tuning the enrichment (and expression) of the metabolisms that are important for the lithification process.

One of the limitations of the taxon‐based prediction of metabolic functions using PICRUSt is that the predictions are biased towards existing reference genomes (Douglas et al., 2020). Nonetheless, as the number of high‐quality reference genomes continues to grow, an increasing number of studies are using, and thus validating, the robustness of this pipeline. Here, it allowed us to predict the metabolic functions of microbialite communities across the five lakes in a rapid and cost‐effective way, providing a general understanding of the metabolisms and associated bacteria that may be important for microbialite accretion.

Importantly, within the context of microbialite formation, this study reinforced the view that ‘who’ is in the community is not as critical as the metabolic potential of the community members and their capacity to respond to biotic and abiotic factors to finely tune the balance between the promotion and dissolution of carbonate precipitation, the net of which ultimately determines the formation of microbialites.

CONFLICT OF INTEREST

The authors declare that there is no conflict of interest.

Fig S1

Click here for additional data file.

Table S1‐S12

Click here for additional data file.

FileS1

Click here for additional data file.

FileS2

Click here for additional data file.

FileS3

Click here for additional data file.

FileS4

Click here for additional data file.