Dynamics of archaea at fine spatial scales in Shark Bay mat microbiomes

The role of archaea in microbial mats is poorly understood. Delineating the spatial distribution of archaea with mat depth will enable resolution of putative niches in these systems. In the present study, high throughput amplicon sequencing was undertaken in conjunction with analysis of key biogeochemical properties of two mats (smooth and pustular) from Shark Bay, Australia. One-way analysis of similarity tests indicated the archaeal community structures of smooth and pustular mats were significantly different (global R = 1, p = 0.1%). Smooth mats possessed higher archaeal diversity, dominated by Parvarchaeota. The methanogenic community in smooth mats was dominated by hydrogenotrophic Methanomicrobiales, as well as methylotrophic Methanosarcinales, Methanococcales, Methanobacteriales and Methanomassiliicoccaceae. Pustular mats were enriched with Halobacteria and Parvarchaeota. Key metabolisms (bacterial and archaeal) were measured, and the rates of oxygen production/consumption and sulfate reduction were up to four times higher in smooth than in pustular mats. Methane production peaked in the oxic layers and was up to seven-fold higher in smooth than pustular mats. The finding of an abundance of anaerobic methanogens enriched at the surface where oxygen levels were highest, coupled with peak methane production in the oxic zone, suggests putative surface anoxic niches in these microbial mats.

Microbial mats in extreme environments provide insights into how physiochemical stresses shape adaptive response in both individual organisms and entire communities. Microbial mats are found globally in environments such as acidic mines drainages 1 , remnant asbestos mines 2 , cold freshwater lakes 3,4 and hypersaline environments in Guerrero Negro 5,6 , Kiritimati Atoll 7 , the Bahamas [8][9][10][11] , and Shark Bay, Australia [12][13][14][15] . Microbial mats in Shark Bay are one of the best examples of modern living microbial mats, and have been proposed to be analogous to Precambrian stromatolites [16][17][18][19] . The extreme conditions of hypersalinity, desiccation stress, and high UV irradiance observed at Shark Bay, may be similar to conditions present on the early Earth [20][21][22][23] . Therefore, it is proposed that living microbial mats can provide us with insights into the evolution of the early biosphere 24,25 .
Microbial mats consist of diverse microbial ecosystems that are formed by multidimensional interactions between biotic and abiotic factors that select for unique functional groups 26 . As a result of these interactions, critical nutrient and biogeochemical cycling (carbon, oxygen, nitrogen and sulfur) usually occur at a millimetre scale and fluctuate during a diel cycle 5,6,27,28 , thus imposing steep chemical gradients within fine spatial resolution 29,30 . This favours niche differentiation among microbial communities, and can also reflect the spatial nutrient requirements and cycling within the mat systems 6 , and potentially how microorganisms interact with each other and the surrounding environment 31 .
Resolving the taxonomic and functional diversity, complexity, and assembly of microbial mat systems has been made possible by employing a combination of classical geomicrobiological techniques with modern next generation sequencing approaches 13,28,32 . However, despite an increasing number of studies indicating archaeal prevalence in some microbial mat systems 6,13,25,33 , the overall role(s) of archaea in microbial mats remains poorly understood. Interestingly, archaea appear to be more abundant in the Shark Bay microbial mats compared to those in open marine ecosystems in the Bahamas 13 . Although early work in Shark Bay indicated microbial mats were dominated by Halobacteria 12,14,34,35 , the techniques used in these studies likely under-represent the true archaeal diversity in Shark Bay mats. Furthermore, new archaeal classifications have been proposed 36,37 , and revision of the classification with closer link to physiologies would help linking functional potential to taxonomic markers such as the 16S rRNA gene, thereby refining and establishing ecological roles of archaea in these environments. Finally, a caveat of many studies relying solely on 16S rDNA data is that there are limits to how much function can be inferred from phylogeny, and thus conducting defined metabolic measurements to support molecular analyses is critical to facilitate a more comprehensive understanding of mat ecosystems.
A recent study deciphered for the first time bacterial community distribution, abundance, and putative interactions at the millimetre scale in Shark Bay microbial mats 28 . Employing iTag deep amplicon sequencing coupled with biogeochemical measurements, the spatial distribution of bacterial groups with mat depth and patterns of niche differentiation were inferred. However, complementary analyses resolving the archaea at a discrete millimetre scale are needed to fully characterise the microbial communities in Shark Bay. Therefore, the aim of this study was to characterise the spatial distribution of archaea with depth combined with defined metabolic measurements, in order to elucidate putative ecological niches and microbial functional roles of this domain in Shark Bay microbial mats. In the current study, deep iTag amplicon sequencing of archaeal small subunit RNA, coupled with biogeochemical measurements (O 2 production and consumption, sulfate reduction, and methane metabolism), were undertaken in microbial mats in Shark Bay. This study has greatly enhanced our understanding of archaeal community dynamics in microbial mat systems, and their putative roles in potentially 'filling the niches' .

Results
Biogeochemical measurements in Shark Bay microbial mats. Smooth mats displayed a shallow O 2 maximum of 321 ± 26 μ M at 1.25 mm depth (average ± standard deviation; n = 3) and a relatively shallow oxycline of 3.2 ± 0.1 mm depth even during peak photosynthesis (light intensity 1780-2040 μ E.m −2 .s −1 ; Fig. 1A). Oxygen penetrated deeper in pustular mats (i.e., down to 6 mm) than in smooth mats, but the maximum concentration of O 2 was two-fold lower (175 ± 14 μ M), with a broad O 2 maximum between 1.5 and 3.75 mm depth (Fig. 1B). The oxygen production and consumption rates were higher in smooth mats compared to pustular and concentrated near the surface (upper 3.5 mm). Similarly, sulfate reduction rates, as evidenced by pixel densities indicating Ag 35 S distribution (Fig. 1D,F), were higher in smooth mats than in pustular mats. Furthermore, the bulk of the microbial activity was concentrated in the upper ca. 3 mm in smooth mats compared to the upper 5-6 mm in pustular mats. Compared to smooth mats, the sulfide concentration at 8 mm depth was 3-7 times lower in pustular mats. The independent metabolic profiles obtained on the mats in the present study corroborated very well with recent measurements taken during the analysis of bacterial populations in these ecosystems 28 . Estimation of richness and diversity. A total of 13,547,552 unfiltered sequences were obtained. To account for sampling bias, databases were subsequently subsampled so that each sample contained 31,000 sequences, with a total of 1,147,000 sequences in all the samples. Rarefaction analysis and non-parametric estimators were performed to determine archaeal diversity and richness. Alpha diversity analyses of all samples were performed at 97% evolutionary distance level of OTUs. Rarefaction analysis showed that the curves have not yet reached an asymptote (Supplementary Figure 1), indicating further sampling is needed to reach saturation of richness. Shannon indices ranged from 4.51 to 5.67 in pustular mats and 3.98 to 5.67 in smooth mats (Supplementary Table 1 Table 2). One-way analysis of similarity test (ANOSIM) indicated the archaeal community structures of smooth and pustular mats were significantly different (global R = 1, p = 0.1%). In particular, archaeal classes MBGB, Thermoplasmata, and Methanomicrobia were abundant in smooth mats but virtually absent (< 0.5%) in pustular mats. The recently proposed Parvarchaeota dominated smooth mats and was also abundant in pustular mats.

Vertical distribution of archaeal communities in smooth mats.
Smooth mat surfaces were highly enriched with Parvarchaeota order YLA114 (42.39% at the surface), followed by Thermoplasmata order E2 (26.37% at the surface) (Supplementary Figure 2). E2 was characterised by archaeal family CCA47, DHVEG-1 (Deep sea hydrothermal vent group I) and Methanomassiliicoccaceae. YLA114 and CCA47 gradually deplete with depth. In contrast, Parvarchaeota order WCHD3-30 and Thermoplasmata lineage DHVEG-1 increase with depth and peak at the bottom layer (15.92% and 24.41% respectively) (Supplementary Figure 2). Crenarchaeota/ Thaumarchaetoa class MBGB peaks at the subsurface (19.42%) (2-4 mm), then gradually decreases with depth. Euryarchaeota was primarily represented by members of class Thermoplasmata, which was also the dominant member at the bottom layer. Hydrogenotrophic methanogens Methanomicrobiales comprise over half of the total methanogenic community, which peaks at the subsurface (2-4 mm) and gradually decreases with depth. Conversely, other known methanogens such as Methanosarcinales, Methanococcales, Methanobacteriales and the recently classified Methanomassiliicoccaceae are methylotrophic, and gradually increase in abundance with depth in smooth mats.
To examine phylogenetic stratification by spatial layers in the mats, cluster and PCA analysis were undertaken. Smooth mats were phylogenetically stratified into distinct spatial mat groups by layer in the mats, primarily surface (0-2 mm) and bottom layers (18-20 mm), and some clustering of subsurface (2-4 mm) and anoxic zone (6-18 mm) (Fig. 2D). STAMP analyses were employed to identify which archaeal members were overrepresented in the mat layers 38 . The surface of smooth mats was enriched in Parvarchaeota order YLA114 (42.39%) Vertical distribution of archaeal communities in pustular mats. In contrast, to smooth mats, pustular mats had > 60% of the archaea represented by Halobacteria, followed by Parvarchaeota ( Fig. 2A). When compared to smooth mats the pustular mats had very low levels of methanogens by sequence analysis (< 0.5%). Halophilic archaeal genera dominanted pustular mats namely, Halogranum (10.48%), Natronomonas (17.41%), Halorhabdus (15.98%), Halobacteriaceae (5.98%) and Halorientalis (0.72%). A candidate genus ArcF12, which belongs to the family Halobacteriaceae, was also found abundant in pustular mats (16.61%) (Supplementary Figure 6). The genus Nitrosopumilus was the dominant Thaumarchaeota in pustular mats, which is an active nitrifier. Pustular mats revealed only two groups in regards to phylogenetic clustering, namely the top layer (0-2 mm) and the subsequent layers (2-10 mm) (Fig. 2B). STAMP analysis revealed that none of the archaea in pustular mats was significantly enriched in any layer.
Methane production. Time series of CH 4 production (nmol CH 4 .cm −3 mat) in smooth mat slurries showed no or only a short (i.e., less than 12 h) lag phase before the headspace concentration of CH 4 started to increase (Fig. 1G,H). Slurries prepared from the oxic part of smooth mats showed the shortest lag phase and accumulated the largest amounts of methane. In contrast, slurries prepared from pustular mats did not produce methane until 12-24 h, suggestive of a smaller active methanogenic population present in pustular compared to smooth mats. When methane concentrations levelled off, potential methane production was calculated from the increase in CH 4 with time. In unamended slurries, the rate of methane production (nmol CH 4 .cm −3 mat.h −1 ) was six to seven times higher in smooth than pustular mats, which correlates with the putative lower abundance of methanogens detected at the molecular level in pustular mats described earlier. The rates were comparable in the oxic and anoxic parts of the smooth mat and only slightly higher in the slurries of the anoxic part of pustular mats compared to those in the oxic part (Fig. 1). In slurries of smooth mats, addition of methanogenic substrates resulted in an increase in methane production rate (nmol CH 4 .cm −3 mat.h −1 ), with the greatest effect of H 2 /CO 2 (Fig. 3). Addition of excess H 2 /CO 2 to smooth mat slurries doubled the rates of CH 4 production and was more effective than addition of excess trimethylamine -a non-competitive substrate.
Inhibition of sulfate reducing bacteria with molybdate increased methane production rates by almost half and inhibition of methanogens using bromoethane sulfonate (BES) or autoclaving reduced activity completely (data not shown). The methane production potential in slurries prepared from the oxic part of smooth mats appeared slightly higher than those in slurries prepared from the anoxic part of these mats (Fig. 1). In pustular mat slurries, the patterns of potential methane production were similar: addition of exogenous substrate increased the rate of CH 4 production up to three times, and inhibition of sulfate reducing bacteria doubled the CH 4 production rate. The rate measured in slurries prepared from the anoxic part of the pustular mats was always higher than in slurries prepared from the oxic part of these mats. When the headspace over the slurries was replaced by air daily, the CH 4 production rate decreased to ca. half in slurries from both oxic and anoxic parts of the smooth mat. In pustular mats, this pattern was the opposite with anoxic slurries having a higher potential methane production rate than oxic slurries.

Discussion
This study describes for the first time the spatial distribution of archaea at high resolution with depth in distinct Shark Bay microbial mats, through a combination of high-throughput 16S SSU sequencing of community DNA and in situ measurements of major metabolic pathways, including oxygenic production and consumption, sulfate reduction, and methanogenesis.

Physiochemical conditions affecting microbial community structures. Smooth and pustular mats
in Shark Bay are found in the littoral zone and are exposed to the atmosphere each tidal cycle. During this period, mats increase in temperature, salinity of the water film overlying the mats increases rapidly, and in certain cases (e.g., high wind), desiccation may set in. Therefore physicochemical conditions inside mats are likely subject to extreme changes at least twice daily. As such, this represents a stressed environment whereby pustular mats are more exposed than smooth mats due the relative location along the tidal transect (Fig. 1). The difference in physicochemical conditions in mat microenvironments -in particular salinity possibly shaping community dynamics -may contribute to the different archaeal diversities observed, as well as the distinct metabolic profiles present in a given mat system (i.e., higher sulfate reduction rates, and higher organic carbon and methane production in smooth mats). The lower overall diversity (albeit higher Halobacterial diversity) observed in pustular mats compared to smooth mats, is thus suggested to be the result of differences in location between the two mat types in the tidal transect of Shark Bay and the use of potentially unique osmoregulation strategies 39,40 . Archaeal groups 'filling the niche' in Shark Bay microbial mats? Previous work on the Shark Bay mats suggested that the Thaumarchaeota may be 'filling the niche' in terms of nitrification 28 , as common bacterial nitrifiers and their pathways are in low abundance in the Shark Bay systems 13,28 . This may be particularly relevant to nitrogen cycling in Shark Bay, as high salinity has been shown to inhibit bacterial nitrifiers in other systems 41 . The detection of the known archaeal nitrifier Nitrosopumilus in the present study may help support this hypothesis, however further investigation comparing archaeal and bacterial nitrification in the Shark Bay systems is needed to definitely delineate the roles of these two groups in nitrogen cycling. Without ammonia oxidisers, accumulation of ammonia is detrimental to microbial ecosystems and therefore Thaumarchaeota such as Nitrosopumilus are proposed to play a role in oxidising ammonia. Potentially this could be facilitated through degradation of glycine betaine 42,43 -an osmoprotectant common in Shark Bay -or deamination of proteins from decaying microbial mat biomass 44 . Whether any microorganisms that possess comammox -a recently discovered pathway that oxidises ammonia to nitrate directly 45 -are present in Shark Bay mats remains to be determined.

Hydrogenotrophic methanogens.
Interestingly, under the conditions of the present study, a wide diversity of methanogen sequences was detected in smooth mats, including orders Methanomicrobiales, Methanosarcinales, Methanobacteriales and Methanococcales. Methanomicrobiales is represented here by the genus Methanofollis, which is strictly hydrogenotrophic and anaerobic 46 . Methanomicrobiales comprises over half of the methanogenic population in smooth mats, suggesting that hydrogenotrophic methanogenesis may dominate methane production in smooth mats. Addition of excess H 2 /CO 2 or trimethylamine to mat microcosms in the present study increased rates of CH 4 production, corroborating the importance of hydrogenotrophic and methylotrophic methane production. The observation of relatively low abundance of methanogenic sequences in pustular mats was supported by the observation that these mats only started producing methane after a lag phase, and that methane production was significantly lower than in smooth mats (Fig. 1). However, inhibition of sulfate reducers stimulated their activity, and thus methanogenesis still likely has a role in pustular mat systems.
The apparent dominance of hydrogenotrophic methanogens in Shark Bay smooth mats is in contrast to findings in other mat systems, where methylotrophic methanogens generally predominate [47][48][49][50] . However, our study agrees with another study that demonstrated methanogenic communities are vertically compartmentalised in hypersaline microbial mats 33 . Methanogenesis is typically dominated by methylotrophic methanogens rather than hydrogenotrophs in hypersaline settings, as sulfate-reducing activities usually outcompete the latter for hydrogen thermodynamically 33,49,[51][52][53][54][55][56] , especially at salinities below 195 ppt 57 . Indeed high sulfate reducing activity was recorded at the surface of Shark Bay smooth mats here and in other studies 28,58 , which one would expect to restrict the hydrogenotrophic Methanomicrobiales. It is hypothesised that the subsurface of Shark Bay smooth mats holds a putative niche -or trophic strategy -for hydrogenotrophic methanogenesis, where Methanomicrobiales peaks potentially due to its proximity to a putative hydrogen source.

Methanogenesis under hypersalinity. Previous studies in Shark Bay mats have described the impor-
tance of the accumulation of osmoprotective molecules, in particular glycine betaine, as adaptive mechanisms to counter the high salinity conditions 13,39,40 . Osmoregulation is likely a critical aspect of Shark Bay mat community ecophysiology. The addition of trimethylamine, a non-competitive substrate for methanogenesis 52,55 and a degradation product of glycine betaine 42 , stimulated methane production 1.6-1.8 fold in smooth mats and 2.4-3.0 fold in pustular mats in the present study. Furthermore, some methanogens in other ecosystems have been shown to even utilise glycine betaine as another alternate non-competitive substrate for methanogenesis 59 . Thus it is hypothesised here that as a result of salinity stress Shark Bay mat communities may need to produce high levels of glycine betaine and other methylated compatible solutes as osmoprotectants, which can then also be potentially utilised for methanogenesis. Although further work is needed to clarify this, such a model would facilitate a competitive advantage for archaea to undergo methylotrophic methanogenesis in these systems.

Co-existence of sulfate-reducing bacteria and methanogens. Despite sulfate-reducing bacteria
(SRB) usually outcompeting methanogens thermodynamically, both co-exist at the subsurface in microbial mats 27,28,57,60 . SRB is represented by the order Desulfobacterales in Shark Bay mats 28 , which was found to be a primary hydrogenotroph in a hypersaline microbial mat in Elkhorn Slough, USA 61 . This group of SRB also peaked at the subsurface (2-4 mm) in smooth mats 28 -and correlated here with measurements of sulfate reducing activity (Fig. 1) -and an early study in Shark Bay also identified that SRB in smooth mats were quantitatively important H 2 -consumers 62 . Though it appears that Desulfobacterales and Methanomicrobiales are competing for H 2 , the SRB and the hydrogenotrophic methanogens might be in a syntrophic relationship. High levels of H 2 are proposed to cause end product inhibition of the carbon sources that methanogens need for growth 33 . Thus, it is proposed here that putative H 2 -consuming Desulfobacterales may act as a regulator to prevent excess H 2 accumulation, enabling the dominance of hydrogenotrophic Methanomicrobiales in smooth mats.
The fluctuation of salinity facilitated by the tidal regime in Shark Bay may also play a role in the co-existence of SRB and methanogens. Salinity measurements taken at different tidal points in Shark Bay illustrate the increase of salinity during low tide (salinity increased beyond 195 ppt) as the mats were subjected to desiccation stress (Supplementary Figure 7). In hypersaline mats of Puerto Rico, methanogens outcompeted sulfate reducers at salinities above 195 ppt 57 . A similar scenario may be prevalent in Shark Bay during low tide, when the high salinity may impede the metabolic activity of SRB, potentially rendering substrate for methanogens. It is proposed here that high salinity may impede the metabolism of SRB, allowing methanogens to take over at the surface, while SRB metabolisms predominate when salinity drops below 195 ppt.
Surface anoxic niche. The high rates of sulfate reduction observed at smooth mat surfaces (Fig. 1) indicate putative microzones of anoxia in oxic parts of the mat, which can further be exploited as possible niches Scientific RepoRts | 7:46160 | DOI: 10.1038/srep46160 for methanogens (shown by high methane production rates at the surface of smooth mats). Conversely, the lower rates of sulfate reduction in pustular mats (Fig. 1), suggested fewer putative anoxic microzones, which was also reflected in the lower methanogen diversity at the taxonomic level, and lower rates of methanogenesis observed. These putative anoxic surface niches were further supported by previous work indicating a high proportion of bacterial anaerobes at the surface of these mats in Shark Bay 28 . Another finding in the present study supporting a potential surface anoxic niche was the detection of abundant Thermoplasmata sequences (order CCA47) at the surface in smooth mats in Shark Bay, and previous studies have primarily detected these groups in oxygen-depleted sediments 63,64 .
These findings corroborate further the suggestion that there are anoxic/suboxic microniches at the surface of smooth mats 28,65 . It has been proposed that tightly regulated metabolisms -potentially mediated by quorum sensing -of microbial consortia such as cyanobacteria, sulfate reducers/sulfur oxidisers, and methanogens/ methanotrophs, may create these oxic and anoxic zones at mat surfaces 66 . Further support for these putative cooperative niches was the finding that phototrophic sulfide-oxidising bacteria (e.g. Chromatiaceae) were significantly enriched at the surface of these Shark Bay mats 28 , in addition to putative methanotrophic bacteria (e.g Verrucomicrobia) in high abundance at mat surfaces 28 . Therefore, we suggest that niche differentiation and metabolic specialisation in the mats may be shaped by microbial interactions between functionally distinct groups (i.e. Methanomicrobiales, Cyanobacteria, Desulfobacterales, Chromatiaceae, Verrucomicrobia) and surrounding physiochemical conditions, leading to the spatially driven unique niches of microbial community structure by layer in Shark Bay microbial mats. The presence of these groups at mat surfaces may be key in efficient cycling of organic C, S, CH 4 , and O in these systems, and the maintenance of these putative anoxic/suboxic surface niches.
Evolutionary consideration. The presence of hydrogenotrophic methanogens described in the present study and high levels of molybdenum (Mo) detected in these modern smooth mats 28 , indicates the systems in Shark Bay might be similar to their ancient counterparts in the Precambrian period 67 . It has been proposed that hydrogenotrophic methanogenesis is an ancient form of methane production 68,69 , whilst enriched Mo was also proposed as an indicator of ancient microbial mats 70 . Furthermore, it has been proposed that the origin of Mo-dependent nitrogenase is linked to a hydrogenotrophic methanogen 71 , suggesting smooth mats may have retained evidence of an ancient origin.
Conclusions. Findings in the present study suggest that archaea may be critical in 'filling the niches' and could be key players in Shark Bay microbial mats. A schematic summarising some of the key interactions inferred from this study is shown in Fig. 4. The present study delineates Parvarchaeota as the dominant archaeal group, though their exact role in ecosystem function is unclear. The smooth and pustular mats have significantly different microbial community structures, and hydrogenotrophic methanogens are enriched despite being potentially less favoured thermodynamically in the presence of sulfate reducing bacteria. Interestingly, it has been suggested that hydrogenotrophic methanogenesis is an ancient form of methane production 68,69 , and thus the modern Shark Bay mats may be good proxies for their ancient counterparts 67 . This high-resolution phylogenetic characterisation, coupled to data describing the spatial distribution of bacteria in these mats 28 , reinforces a putative model of a surface anoxic niche in these systems. It should be acknowledged that there are limitations to the present study, given the focus at the molecular level was on one taxonomic marker (16S rDNA of archaea) from total community DNA. To this end, future functional metagenomic and metatranscriptomic profiling at different depths (also targeting cDNA and active communities), coupled with microbial activity measurements (CH 4 , O 2 ), along with techniques such as Nano FISH-SIMS, stable isotope analyses, and metaproteomics, are required to confirm the data and the putative niches and interactions proposed. In addition, although the differences of metabolic measurements between surface oxic and deeper layers is highly reproducible some of the differences between layers is small, and thus follow up studies over different depths (and seasons) are needed to ascertain true functional roles in these proposed niches. Nevertheless, our data suggests there is likely a reservoir of untapped archaeal diversity in these systems, and metabolic cooperation between key microbial groups is proposed to be important in efficient cycling of key nutrients.

Sample collection. Smooth and pustular mats were sampled on 16 th June 2013 from Nilemah, Hamelin
Pool, Shark Bay (26°27′336″S, 114°05.762″E), using methods as previously described 12,28 . Of the two mats, pustular mats were located closer to the shore than smooth mats (Fig. 1). At the time of sampling, water temperature was 24.8 °C, salinity 68 PSU, and pH 8.13. Triplicate mat samples were taken using a sterile scalpel and placed in sterile containers. Half of the mat samples were immediately preserved in RNALater (Ambion, Life Technologies) for downstream DNA analyses. In order to determine the archaeal composition in successive mat layers, the samples were dissected into 2 mm intervals using sterilised blades as previously described 28 . Measurements of key biogeochemical properties. In order to compare major microbial processes between the two mat types, and to determine the approximate depth of the oxic and anoxic zones, depth profiles of oxygen, sulfide and pH were measured in small mat samples (~5 × 5 cm) submerged in 3 cm of water collected from the site. Profiles were measured in triplicate using needle microelectrodes 27 used in combination with a Unisense PA 2000 picoammeter, and pH and ion-selective sulfide needle (Microscale Measurements, The Netherlands) were deployed with a high-impedance microscale measurements millivolt meter (Microscale Measurements). The overlying water was gently stirred during the measurements. Because of the variations in mat surface topography and slight variations in the thickness of the different mat layers in smooth mats and thickness of the pustules, the depth of oxygen penetration and the observed maximum values of O 2 and sulfide varies. Replicate measurements were made ca. 1 mm apart along a linear transect. This transect was marked and upon completion of electrode measurements, dissected for two-dimensional mapping of sulfate reduction activity (see below). Light measurements were performed using a LiCor LI 250 meter equipped with a SA190A quantum sensor.
Microbial metabolic rates were estimated using the light dark shift (for O 2 production and consumption 72,73 ). Briefly, the O 2 concentration with depth was determined in the light. Following, the mats were incubated in the dark and O 2 profiles measured every 10-15 min for ca. 1.5-2 h until profiles resembled dark profiles measured at 0400 h (i.e., depicting O 2 diffusion into the mat). Consecutively, mats were exposed to natural sunlight and O 2 profiles measured as before until the original light profile was approached. Oxygen consumption and production were calculated from the decrease, and increase minus decrease in [O 2 ], respectively. Oxygenic photosynthesis produces oxygen, and aerobic respiration, chemolithotrophic sulfide oxidation and chemical reactions involving O 2 are the key contributors to oxygen consumption. The 2D-distribution of sulfate-reducing activity was visualized using small strips (ca. 5 × 10 mm) of ref. 35  . Upon reduction by SRB, the metabolic product 35 S 2− binds with the Ag foil surface and the radioactivity can be visualized using a radioactivity gel scanner 74 . The result is a quantitative high-resolution map of the distribution of SRB activity with depth. Freshly collected mat sections were cut vertically and placed on the Ag-foil strips. After 4 and 8 h of incubation at ca. 30 °C for the top and bottom of the core, respectively, the mat samples were removed and the remaining 35 SO 4 2− rinsed off the foil using distilled water. The foils (containing 35 S 2− preserved as Ag 35 S during SR) were kept in the dark and scanned using the BioRad Molecular Imager System GS-525 (Hercules, CA) to visualise a 2-D Ag 35 S distribution. The result is a digital pixel map, in which darker pixels represent microzones of higher sulfate reducing activity. The spatial resolution of this technique is 140 × 140 μ m.
Methane production was measured in samples upon return to the laboratory 56,60 . Based on in situ microelectrode O 2 profiles, mats were separated in an oxic and anoxic part (0-3 and 3-20 mm for smooth mats and 0-6 and 6-10 mm for pustular mats, respectively). Slurries were prepared from entire layers of the mat, containing both organic (e.g., biomass and EPS) and inorganic (e.g., mineral) fractions. Individual slurries may comprise slightly different amounts of biomass, however identical biomass is not needed as the focus of these assays is on the overall differences in total community activity (methane production) between two different mat systems. Biomass differences may contribute to some of the differences observed and this is ecologically relevant. The samples were homogenised in 1:1 (v/v) mixture of 0.22 μ m-filtered site water and mat and transferred to 18-ml crimp seal bottles (12 ml slurry and 6 ml headspace). Manipulations of mats and slurries were performed under a N 2 atmosphere. The headspace of individual vials was flushed with N 2 for 5 min at T = 0 except for the treatment in which the headspace was flushed with air daily. Triplicate samples were incubated in the dark at 25 °C on a shaker table (~60 rpm). The evolution of methane in the slurries was followed for 7-10 days by headspace analyses using gas chromatography with flame-ionisation detection (Shimadzu GC14A) and a PoraPak Q column and PeakSimple integration software 60 . In addition to appropriate controls, the following amendments were made to the slurries: trimethylamine (10 μ M final concentration), H 2 /CO 2 (75/675 μ M), molybdate (28 mM) and BES (5 mM). In one treatment, the headspace was replaced daily by air.
Nucleic acid extraction. Total community genomic DNA was extracted in duplicate from each microbial mat layer (10 layers from smooth mat, 5 from pustular mat) employing the MoBio PowerBiofilm DNA Isolation Kit (MO BIO Laboratories, Carlsbad, USA) according to the manufacturer's instructions. A cross section of the mats indicating where they were layered is shown in Supplementary Figure 8. The concentrations and purity of extracted DNA were determined spectrometrically, and the quality checked by PCR amplification of archaeal 16S rRNA genes 12,13,28 . 16S SSU amplicon sequencing and analysis. Paired-End sequencing of archaeal 16S SSU genes from each mat layer in duplicate was performed using a MiSeq V3 2 × 300 bp kit on an Illumina MiSeq desktop sequencer 76 . Amplicons were produced using barcoded Arch-304F-Nex and Arch-915R-Nex primers. Amplicon primers were Arch-304F-(5′ TCGTCGGCAGCGTCAGATGTGTATAAGAGACAGCCCTAHGGGGYGCASCA-3′ ), Arch-519R-(5′ -GTCTCGTGGGCTCGGAGATGTGTATAAGAGACAGGWGCYCCCCCGYCAATTC-3′ ). MOTHUR software package version 1.33.0 77 was employed. Quality control steps include removing sequence reads longer than 500 bp, removing homopolymers longer than 8 bp, and removing those containing ambiguous nucleotides. Sequences were then aligned (Needleman-Wunsch pairwise alignment) with the Silva bacteria database version 119 78 , and filtered to remove non-informative columns. Subsequently, the UCHIME function included in MOTHUR was employed to remove chimeric sequences in de novo mode 79 . After chimeras were removed, sequences were classified using default settings in MOTHUR against the Greengenes database 80 (August 2013 released version). SSU sequences from each group were then randomly subsampled to 31,000 sequences (the fewest among the dataset) in order to normalise the sequence number for downstream statistical analyses. Subsequently, the sequences were clustered de novo into operational taxonomic units (OTU) at a genetic divergence level of 3% 81 (or at 97% sequence identity). Taxonomic richness and diversity estimators were calculated using MOTHUR 77 .
Statistical analyses. PRIMER 6 plus permanova packages were used for multivariate analyses, in which the relative abundance of each individual taxa was generated. The data was relative abundance of clustered OTUs at 97% cut-off for 16S rDNA based on archaea. The average relative abundance determined between duplicate layers is shown in Supplementary Table 2. Unless otherwise noted, abundance of OTUs were square root or natural log transformed in order to de-emphasize the large values; hence rare sequences of low abundance can also be taken into account. One-way Analysis of similarity test (ANOSIM) was performed to examine whether the microbial communities between both mat types were significantly different. PCA biplot analysis was performed in the Factomine package in R, to indicate how different archaeal groups putatively correlate to different depths in the mats. To further explore microbial stratification, statistical analysis of metagenomic profiles (STAMP) was undertaken 38 . Taxonomic groups from one layer were compared to the same taxonomic groups in all other layers, to identify which archaeal taxa in each layer were best discriminated significantly between other designated layers. Welch's t-test was employed to reduce type I errors 82 , and Benjamini-Hochberg (i.e. false discovery rate) correction was applied as a multiple-hypothesis test correction. All quoted P-values represent corrected values (equating to q), with values < 0.05 considered significant 83 . A filter was applied to remove taxa with a q value > 0.05.