Metatranscriptomic analysis of a high-sulfide aquatic spring reveals insights into sulfur cycling and unexpected aerobic metabolism

Zodletone spring is a sulfide-rich spring in southwestern Oklahoma characterized by shallow, microoxic, light-exposed spring water overlaying anoxic sediments. Previously, culture-independent 16S rRNA gene based diversity surveys have revealed that Zodletone spring source sediments harbor a highly diverse microbial community, with multiple lineages putatively involved in various sulfur-cycling processes. Here, we conducted a metatranscriptomic survey of microbial populations in Zodletone spring source sediments to characterize the relative prevalence and importance of putative phototrophic, chemolithotrophic, and heterotrophic microorganisms in the sulfur cycle, the identity of lineages actively involved in various sulfur cycling processes, and the interaction between sulfur cycling and other geochemical processes at the spring source. Sediment samples at the spring’s source were taken at three different times within a 24-h period for geochemical analyses and RNA sequencing. In depth mining of datasets for sulfur cycling transcripts revealed major sulfur cycling pathways and taxa involved, including an unexpected potential role of Actinobacteria in sulfide oxidation and thiosulfate transformation. Surprisingly, transcripts coding for the cyanobacterial Photosystem II D1 protein, methane monooxygenase, and terminal cytochrome oxidases were encountered, indicating that genes for oxygen production and aerobic modes of metabolism are actively being transcribed, despite below-detectable levels (<1 µM) of oxygen in source sediment. Results highlight transcripts involved in sulfur, methane, and oxygen cycles, propose that oxygenic photosynthesis could support aerobic methane and sulfide oxidation in anoxic sediments exposed to sunlight, and provide a viewpoint of microbial metabolic lifestyles under conditions similar to those seen during late Archaean and Proterozoic eons.


INTRODUCTION
The combination of molecular techniques and high-throughput sequencing in microbial ecology has revolutionized the ability to probe questions regarding microbial community structure and function (Handelsman, 2004). Metagenomic and metatranscriptomic approaches have provided valuable insights into the roles of uncultivated lineages within a microbial community (Pelletier et al., 2008;Liu et al., 2012;Kozubal et al., 2013;Sheik, Jain & Dick, 2014) as well as in delineating the relative importance and contributions of various microbial taxa to observed geochemical processes occurring in habitats characterized by complex biodiversity, such as soil (Urich et al., 2008;Shrestha et al., 2009), surface and deep ocean water (Poretsky et al., 2005;Frias-Lopez et al., 2008;Poretsky et al., 2009;Baker et al., 2013), deep ocean hydrothermal vent systems (Xie et al., 2010;Lesniewski et al., 2012), and hot springs (Liu et al., 2011;Burow et al., 2013). In complex habitats like these, small-subunit (SSU) rRNA gene-based surveys and metagenomic approaches can lead to exhaustive lists of taxa detected and speculation of the potential roles of dominant taxa in metabolic processes. As such, a metatranscriptomic approach better allows for the determination of microbial populations that are actively involved in processes of interest.
In particular, metatranscriptomics represents an ideal tool to study sulfur cycling processes and the microbial lineages involved at Zodletone spring, a sulfide and sulfur-rich spring in southwestern Oklahoma. The spring itself is a sunlight-exposed shallow aquatic habitat with highly reduced sediments, in which O 2 is undetectable (Senko et al., 2004;Bühring et al., 2011). Because of the lack of detectable O 2 , it has also been suggested that photosynthesis by Cyanobacteria in spring's sediment is anoxygenic (Bühring et al., 2011). The water overlying the sediments has low O 2 concentrations (2-4 µM) (Bühring et al., 2011), high levels of methane from the subsurface, and high concentrations of sulfur in various forms (Donovan, Younger & Ditzell, 1988).
Sulfur cycling processes at the site have been described as anaerobic. Anoxygenic sulfide oxidation appears to be light-driven and anaerobic (Elshahed et al., 2003;Senko et al., 2004), and indeed, Chromatiales and Chlorbi (the anaerobic purple and green sulfur bacteria) lineages have been detected among SSU rRNA-based surveys (Elshahed et al., 2003;Youssef, Couger & Elshahed, 2010); however, it is unknown whether these lineages are the sole contributors to sulfide and sulfur oxidation or if other lineages, such as Epsilonproteobacteria, which has been detected in the spring (Elshahed et al., 2003;Youssef, Couger & Elshahed, 2010;Bühring et al., 2011), might also be involved. Similarly, sulfate reduction also appears to be a dominant process in the anoxic sediments during dark periods (Elshahed et al., 2003;Senko et al., 2004), but while many anaerobic Deltaproteobacteria lineages known to be involved in sulfate reduction have been detected (Elshahed et al., 2003;Youssef, Couger & Elshahed, 2010), it is unclear if other sulfur-reducing processes (e.g., sulfur reduction, thiosulfate reduction or disproportionation) also occur in situ and which microbial lineages are involved. Indeed, Zodletone spring sediments harbor an extremely diverse microbial community (Elshahed et al., 2003;Youssef, Couger & Elshahed, 2010), but SSU rRNA-based surveys alone cannot tell us specifically which of the taxonomic lineages detected (especially those with no or few members that have been validly described) are actively involved in any of the sulfur cycling processes, or other photo-and chemotrophic processes, occurring at the site.
Thus, a metatranscriptomic approach, especially if coupled to a detailed geochemical characterization, could be extremely useful. In this study, we used metatranscriptomics to examine the roles of various microbial lineages in processes occurring in Zodletone spring sediments throughout the course of a day, with the goal of identifying microorganisms involved in general photo-and chemotrophic processes as well as sulfur cycling processes. The data presented here highlight taxa involved in a complex sulfur cycle and also provide some unexpected findings regarding the production and utilization of oxygen in the anoxic sediments at the spring's source.
The soft, black, biomass-rich sediments at the spring's source are immobile and represent a stable environment. Although it was previously reported that spring water is anoxic (Senko et al., 2004), a more recent study using microsensors found that though oxygen was not detectable within and directly above (<1 mm) the source sediment, concentrations increased to 2-4 µM at 2 mm above the source sediment-water interface; therefore, source sediments and overlying waters can be classified as anoxic and micro-oxic (<5 µM), respectively (Bühring et al., 2011).

Sampling description
Zodletone source sediment and groundwater were sampled in early November 2009 over a 24 h period to catch the diurnal fluctuations in sulfur cycling processes previously observed, i.e., that sulfide oxidation and production appear to occur during sunlit vs. dark periods, respectively (Elshahed et al., 2003;Senko et al., 2004). Conditions on the day of sampling were sunny and mild, with air temperatures recorded as follows for each time of sampling: 24.4 • C at 12:15 (daytime; peak sunlight), 22.8 • C at 17:15 (early evening; one half hour before sunset), 15.6 • C at 22:15 (night; full darkness), and 15.0 • C at 07:30 (early morning; one half hour after sunrise). At each of these time points, source sediment (approx. 5 g) was aseptically collected, preserving transcripts by immediate (within five seconds) on-site immersion of sediment in 9 mL LifeGuard TM Soil Preservation Solution (MoBio Laboratories, Inc., Carlsbad, CA). Samples were mixed well then frozen by submersion in an ethanol bath super-cooled with dry ice. Using this procedure, each sediment sample-LifeGuard mixture was frozen in under one minute. Each sample was collected from the top layer (5 cm) of sediment from the south-facing corner of the spring (Fig. S1), directly adjacent to the spring's source of water (as identified by visible bubbling), where the sediment has a uniform appearance, characterized as soft, black, with few to no leaves present. Care was taken to ensure minimal disturbance, and subsequent samples were obtained, using the same techniques, from the same location and depth immediately adjacent to the original site sampled. All sediment samples were transported on dry ice to our laboratory and stored at −80 • C for RNA extraction.
Water from the spring's source was also collected for various geochemical analyses. For pH and anion analysis, water from the source was collected, filtered through 0.2 um filter, and frozen immediately using the dry ice/ethanol bath. For sulfide preservation, water samples (2.5 ml) were collected and directly injected into serum tubes containing 2.5 ml anoxic 10% zinc acetate (N 2 headspace). Water samples were transported either on dry ice (tubes for pH/anion analysis) or on ice (tubes for sulfide analysis) to our laboratory and stored at −20 • C (tubes for pH and anion analyses) or 4 • C (tubes for sulfide measurements). Water samples were collected, transported, and stored in a similar matter for pH and sulfide analyses again in August 2014.

Analytical methods
From water samples, pH was measured under N 2 :CO 2 (80:20) headspace using a pH electrode. Anion (chloride, nitrate, nitrite, phosphate, sulfate, sulfite, and thiosulfate) concentrations were determined by ion chromatography (Dionex, model DX500 fitted with the AS-4A column; Dionex Corporation, Sunnyvale, CA). To minimize oxidation of thiosulfate and sulfite, standards and sample dilutions into Dionex vials were prepared in an anaerobic glove bag (Coy Laboratory Products, Inc., Grass Lake, MI). Sulfide was measured by the methylene blue assay (Cline, 1969).

RNA extraction, processing, and pyrosequencing
RNase-free materials and reagents were used in all RNA extraction and processing steps. RNA was extracted from sediment within two months of sampling (RNA PowerSoil ® Total RNA Isolation kit; MoBio). Remaining DNA was then digested (RNase-free DNase-I; MoBio), and an enzymatic rRNA removal step was performed (mRNA-ONLY Prokaryotic mRNA Isolation Kit; EPICENTRE Biotechnologies, Madison, WI). Samples were stored at −80 • C until a second rRNA-removal step, based on subtractive hybridization, was performed (MICROBExpress TM Bacterial mRNA enrichment kit; Ambion, Austin, TX). Enriched mRNA was amplified via in vitro transcription of synthesized cDNA (MessageAmp TM II-Bacteria Prokaryotic RNA Amplification Kit; Ambion). We quantified total RNA (ng/µl) before and after ribodepletion steps. As a control, RNA was also amplified from the non-ribodepleted RNA (RNA prior to rRNA depletion steps) fraction of one sample (12:15 time point) in order to examine the efficiency of rRNA removal, to examine the microbial community structure of the non-ribodepleted sample, based on analysis of SSU rRNA gene transcripts present, and to describe the effects of ribo-depletion on SSU rRNA-based microbial community structure.
Amplified RNA was assessed for quality and quantity by both denaturing gel electrophoresis and spectrophotometric analysis (NanoDrop Products, Wilmington, DE). The quality of the 17:15 was too poor for sequencing; thus, library preparation and sequencing via Roche Titanium 454 pyrosequencing was conducted from the remaining three samples (12:15, 22:15, and 07:30), using the service of a commercial provider (MOgene, LC., St. Louis, MO).

Sequence handling and analyses
(i) Sequence data handling. Raw sequence data (.fna and .qual files) were uploaded to MG-RAST (v2) (Meyer et al., 2008) for data handling and functional annotation. From these total reads datasets, several data subsets were also created. Microbial SSU rRNA data subsets were created by downloading reads from MG-RAST with significant alignment (min. alignment length = 50 bp; max. E-value = 1E-30) to sequences from within the Greengenes (Bacteria) and SILVA SSU (Archaea) databases. Potential mRNA data subsets from each ribodepleted sample (12:15, 22:15, and 07:30) were created by removing all sequences with significant similarity (max. E-value = 10E-3) to bacterial, archaeal, or eukaryotic large and small subunit ribosomal RNA sequences, as identified using the BLASTN alignment tool against the nt database (Altschul et al., 1997). As well, reads from potential mRNA data subsets were assembled for more detailed analysis of sulfur and oxygen cycling genes. To assemble the raw reads, we performed three separate assemblies using Newbler, the 454 assembly software, on untrimmed reads, reads trimmed at 84, and reads trimmed at 63 cycles. The three different trimming lengths were used to reduce the number of artificial contigs produced due to poor qualities at the end of the contigs (Wiley et al., 2009). The contigs from the three separate assemblies were combined and further assembled with Phrap (Ewing & Green, 1998). Potential mRNA and assembled datasets were also uploaded to MG-RAST (v2) for analyses.
(ii) Taxonomic classification of SSU rRNA reads. Taxonomic assignment of each SSU rRNA read was performed using the on-line Classifier tool available from the Ribosomal Database Project (Wang et al., 2007). In order to characterize the effect of the rRNA removal steps on microbial community analysis, microbial community composition (characterized by the number of taxa detected at the phylum, class, and order levels) of the 12:15 total RNA control sample was compared to that of its (12:15 time point) corresponding ribodepleted dataset. Population distribution at each taxonomic level was also determined using relative frequencies of each taxon, and the two samples were statistically compared by a chi-square test based on an r × k contingency table.
(iii) Taxonomic assignment of mRNA transcripts. A BLASTX-MEGAN based analysis was used to assign NCBI taxonomic affiliations to mRNA transcripts. First, sequences from potential mRNA data subsets were examined using the BLASTX search query against the NCBI nr database; results (top 5-10 most closely related hits with a max. E-value cutoff of 10E-5) were saved and summarized in searchable tab-delimited text files. The BLASTX results produced were then used by MEGAN, v4 (Huson et al., 2011) to assign a taxonomic affiliation to each mRNA transcript. Sequences in which top BLASTX hits identified with different taxonomic lineages (i.e., sequences with similar functions but belonging to different phyla) typically remained "taxonomically unclassified" by the program, thus minimizing the chances of false classifications that might occur if only the top BLASTX hit had been used.
(iv) Functional analyses of mRNA transcripts. Potential mRNA data subsets were analyzed in MG-RAST, v2, for functional gene annotation and SEED subsystem determination (Meyer et al., 2008), using the following parameters: minimum alignment length, 34 bp; maximum E-value, 0.01. A number of rRNA sequences remained in these datasets; thus, the final number of potential mRNA reads for each dataset was determined by subtracting the number of sequences with significant alignments (min. alignment length = 50 bp; max. E-value = 1E-30) to SSU and large subunit (LSU) rRNA databases from the total number of sequences. Functional gene assignments from each subsystem were downloaded as sortable, searchable spreadsheets.
Additionally, individual reads and assembled contigs coding for genes putatively involved in sulfur and oxygen transformations were identified by mining all data files (MG-RAST annotation spreadsheets, MEGAN files and BLASTX summary files) for genes identified in MetaCyc pathways (Caspi et al., 2008) or by primary literature as being involved in metabolisms of interest, and subjected to further in-depth analysis.

Nucleotide sequences accession numbers
Sequences for each of four samples were deposited to GenBank through the sequence read archive and can be retrieved from the following accession numbers: SRX025760, SRX025825, SRX025796 [corresponding to ribodepleted 22:15 (night), 07:30 (early morning), and 12:15 (afternoon) samples, respectively], and SRX331735 (corresponding to the 12:15 non-ribodepleted RNA control sample). Various datasets are publicly available in MG-RAST (v3)

Geochemistry of Zodletone spring water
Sulfide measurements from source water samples at the time of sampling (Nov, 2009) and again in Aug, 2014, indicated a progressive decrease in sulfide values during the day, followed by increase in sulfide levels overnight (Fig. 1). These results confirm previous findings from two previous studies that utilized either direct geochemical measurements at the spring source or laboratory incubations of spring source sediments, and demonstrated that sulfate production from sulfide is primarily a light-dependent process, whereas sulfate reduction/sulfidogenesis occurs predominately nocturnally/in dark incubations (Elshahed et al., 2003;Senko et al., 2004). The results of these studies conducted over more than a decade, along with the documentation of a similar pattern in this study (Fig. 1), demonstrate that the observed sulfur cycling dynamics in Zodletone spring source water is an ecologically relevant and temporally stable process, with sulfur oxidative processes occurring during the day and reductive processes increasing at night.

Pyrosequencing dataset information
Pyrosequencing of four samples (ribodepleted 22:15, 07:30, and 12:15 samples, and the 12:15 non-ribodepleted 12: 15 RNA sample) yielded a total of 1,122,073 reads, with average read lengths ranging between 344-427 nucleotides (Table 1). Ribodepletion steps resulted in a two-to three-fold enrichment in the percentage of potential protein-coding, or mRNA reads, which comprised 11.0-17.7% of total RNA reads in the ribodepleted datasets (Table 1). From potential mRNA subsets, a total of 30,887 reads (22.2-26.6% of total potential mRNA reads) were functionally annotated by MG-RAST (Table 1). BLASTX analysis of these datasets, however, resulted in 44.7-48.9% of potential mRNA transcripts having significant (maximum E value = 1E-5) hits against the nr database (Table 1). Finally, 10.2-13.4% of the potential mRNA transcripts were taxonomically unclassified (could not be placed into any domain or phylum with certainty) (Table 1).

Taxonomic classification of SSU rRNA reads and the effect of rRNA depletion on microbial community structure
Microbial community structure from each SSU rRNA data subset was evaluated by phylogenetic analysis of SSU rRNA reads. Results from our non-ribodepleted sample showed a remarkably similar microbial community structure to what has been described in the past (Elshahed et al., 2003;Elshahed et al., 2004;Youssef, Couger & Elshahed, 2010), with dominant phyla including Proteobacteria, Bacteroidetes, Firmicutes, Actinobacteria, , and sulfate (SO 2− 4 ) concentrations, respectively. Y-error bars represent the standard deviation of duplicate measurements for sulfate and triplicate measurements for thiosulfate. Sulfite, pH, and sulfide values were based on single measurements only. Source water was sampled again in August 2014 to measure pH and sulfide in triplicate and verify patterns of sulfide loss and generation during the day and night, respectively (B). Spirochaetes, Euryarchaeota, and Chloroflexi. Candidate divisions (e.g., BRC1, OP11, OD1, WS3) and several "rare phyla" (i.e., taxa that comprise <0.1% of total microbial community), defined in a previous study from Zodletone (Youssef, Couger & Elshahed, 2010), were also represented in SSU rRNA subsets ( Fig. S2 and Table S2). For example, although the phylum Planctomycetes had been previously determined to comprise 0.01-0.1% of the microbial community in Zodletone source sediment (Youssef, Couger & Elshahed, 2010), Planctomycetes-affiliated SSU rRNA sequences comprised 0.563% of the total SSU rRNA sequences from our non-ribodepleted SSU rRNA dataset and 1.45-2.15% of sequences from ribodepleted datasets SSU rRNA (Table S2). The removal of ribosomal RNAs did not appear to reduce the number of taxa detected at phylum, class, or order levels (Tables S1 and S2). For example, 31 phyla were detected in both ribodepleted and non-ribodepleted RNA samples from 12:15 (Tables S1 and S2). However, though similar taxa are detected from non-ribodepleted vs. ribodepleted samples (Fig. S2), the two microbial community structures, based on rel. frequencies of each taxon, were significantly different at the phylum, class and order levels (p < 0.001). Thus, as the ribodepletion steps inherently increase the error associated with comparing the relative abundance of one taxonomic group to that of others within a sample, ribodepleted datasets should be limited to only the detection of active microbial taxa in a sample and not to estimate taxon abundance or to compare microbial community structures of two or more different samples. Therefore, the results presented hereafter in this manuscript will focus on mRNA transcripts from our ribodepleted samples.
With respect to the phototrophic lineages, Cyanobacteria-affiliated mRNA transcripts were most abundant, especially in the early morning (07:30) sample, comprising 1.26% of total potential mRNA transcripts at that time ( Fig. 2 and Table S3). In the afternoon (12:15), transcripts affiliated with the class Chloroflexi were most abundant, comprising 0.637% of total potential mRNA transcripts (Table S3). Protein-coding transcripts mapping to Chlorobi (green sulfur bacteria) were also detected, but were less abundant (avg abundance = 0.181%) than those mapping to Cyanobacteria or Chloroflexi (Table S3). Protein-coding transcripts that mapped to the phototrophic purple sulfur bacteria (Chromatiales) and the purple non-sulfur bacteria (Rhizobiales, Rhodospirillales, Rhodobacterales, and Rhodocyclales) were also detected among potential mRNA datasets at all three time points (Table S3).
With respect to lineages typically composed of aerobic chemotrophs, mRNA transcripts mapping to the predominantly aerobic chemoheterotrophic lineages Actinobacteria and Halobacteria were more highly transcribed at 12:15 (comprising 2.77 and 0.497% of total potential mRNA transcripts, respectively) compared to the other two time points (Fig. 2 and Table S3). Similarly, mRNA transcripts assigned to lineages containing aerobic sulfide/sulfur-oxidizing chemolithotrophs, such as the Sulfurovum group within the order Campylobacterales in the class Epsilonproteobacteria (Inagaki et al., 2004), and the order Thiotrichales within the Gammaproteobacteria (Garrity, Bell & Lilburn, 2005), were also detected at their highest relative abundance values at 12:15 ( Fig. 2 and Table S3). Protein-coding transcripts mapping to aerobic type I methanotrophs, i.e., Methylococcales within the Gammaproteobaceria (Bowman, 2005), however, were approximately 5-fold more abundant at night (22:15) compared to 07:30 and 12:15 datasets (Table S3).
Protein-coding transcripts assigned to taxa typically composed of anaerobic chemotrophs, including the Bacteroidetes, Spirochaetes, Fibrobacteres, Deltaproteobacteria, and Methanomicrobia, were more highly represented in the early morning (07:30) dataset compared to the other two time points (Fig. 2 and Table S3). Protein-coding transcripts assigned to classes Clostridia and Bacilli within Firmicutes were most abundant at 12:15 ( Fig. 2 and Table S3).

Functional distribution of mRNA transcripts
Based on MG-RAST analysis, carbohydrate metabolism and protein metabolism were among most abundant subsystem categories detected (Fig. 3A). Specifically, xylose utilization was the most abundant among metabolic subsystems detected (Fig. 3B); utilization of xylose and other simple sugars (e.g., arabinose, maltose, and ribose) were more abundant at 07:30 and 12:15 than at 22:15 (Fig. 3B). Subsystems related to photosynthesis, i.e., Photosystems II and I, and light-harvesting protein complexes known as phycobilisomes (De Mot et al., 1999), were most highly transcribed in the early morning (07:30) (Figs. 3A and 3B), consistent with the abundance of Cyanobacteria-affiliated mRNA transcripts at that time ( Fig. 2 and Table S3).
Transcripts coding for the copper-containing particulate methane monooxygenase, one of the two major enzymes involved in methane oxidation (Nguyen et al., 1998;Schwarz, 2001), were highly abundant at 22:15 (Fig. 4), consistent with the high relative abundance of Methylococcales-affiliated mRNA transcripts also observed at this time (Table S3). Methanogenesis transcripts, one the other hand, were more highly abundant in the early morning (07:30) sample (Fig. 3B), which is consistent with the higher relative abundance of Methanomicrobia-affiliated mRNA transcripts also detected at that time (Table S3). Genes for hydrogenases (Semrau et al., 1995) and cellulosomes, which compartmentalize cellulose degradation in anaerobes (Vignais & Billoud, 2007), were also most abundant at 07:30 (Fig. 3B).
Bacterial proteasomes, which compartmentalize protein degradation enzymes and processes inside a bacterial cell (Westley & Green, 1959), along with protein degradation subsystems were detected at a higher abundance at 12:15 compared to the other two time points (Fig. 3B). Genes for fermentation pathways (butyrate and lactate) and sulfur oxidation, also appeared to be more highly transcribed at 12:15 compared to the other two timepoints (Fig. 3B).

Diversity and abundance of S cycling transcripts
(i) Oxidation of sulfur compounds. Consistent with sulfide loss occurring during the day (Fig. 1), transcripts coding for sulfide oxidation were more abundant (0.034% of potential mRNA reads, respectively) at 12:15 (Fig. 5) compared to other time periods. Surprisingly, many transcripts coding for sulfide oxidation (flavocytochrome c sulfide dehydrogenase) mapped to the phylum Actinobacteria (Fig. 5). A drop in thiosulfate was also observed during the day (Fig. 1); likewise, thiosulfate oxidation transcripts belonging to the Sox system were also observed at their highest abundance (0.056% of potential mRNA reads) at 12:15 (Fig. 5). A closer look at Sox genes transcribed at 12:15 showed that the all mapped to Gammaproteobacteria (Fig. 5); specifically, mRNA transcripts assigned as Thiotrichales, which includes chemolithotrophic sulfur-oxidizing bacteria (Garrity, Bell & Lilburn, 2005), comprised 65.0% of these, and those assigned as Chromatiales, which are the anoxygenic sulfide-oxidizing phototrophs known as the purple sulfur bacteria (Imhoff, 2005), comprised 25.0%. Gene transcripts involved in sulfur/polysulfide-oxidation processes were also more abundant at 12:15 (0.022% of potential mRNA reads) than at other time points; the vast majority of these genes expressed at all time points (91.7%) remained taxonomically unclassified (Fig. 5). Chlorobi (green sulfur bacteria)-affiliated transcripts only comprised a small number of sulfide and polysulfide oxidizing transcripts (Fig. 5).
(ii) Reduction of sulfur compounds. Consistent with patterns of sulfide production occurring overnight (Fig. 1), transcripts coding for many genes involved in pathways involving the reduction of sulfur compounds were most abundant at night (22:15) and in the early morning (07:30) (Fig. 5). Deltaproteobacteria-and Firmicutes-affiliated transcripts comprised a large proportion of dissimilatory sulfate-reduction genes, which were most highly transcribed at night (0.034% of potential mRNA reads) compared to other time points. Deltaproteobacteria-affilated transcripts also comprised a large majority of genes transcribed for sulfur compound reduction, which were most highly transcribed at 07:30 (0.019% of potential mRNA reads) (Fig. 5). Transcripts coding for tetrathionate reduction were also detected at their highest levels at 07:30 (0.015% of potential mRNA reads); of these, many also mapped to the Deltaproteobacteria (Fig. 5). Interestingly, the majority of transcripts involved in dissimilatory sulfate reduction and tetrathionate reduction detected at 12:15 were taxonomically unassigned compared to those detected either at 22:15 or 07:30 (Table 2 and Fig. 5).

Expression of genes coding for O 2 production and utilization processes in source sediment
Mining our potential mRNA data subsets for reads and assembled contigs involved in oxygen cycling revealed the diversity of microbes participating in both oxygen production and utilization in Zodletone source sediment (Fig. 4). Cyanobacteria, which were previously considered anoxygenic in the source (Bühring et al., 2011), showed active transcription of the psbA gene (Fig. 4), which codes for the oxygen-producing P680 reaction center D1 protein in Photosystem II (Zouni et al., 2001;Kamiya & Shen, 2003;Ferreira et al., 2004). In addition to cyanobacterial genes, eukaryotic psbA genes (classified as "unclassified Eukaryota" in Fig. 4) had top BlastX hits to psbA genes found among different diatom genomes, including Thalassiosira pseudonana, Odontella sinensis, and Phaeodactylum tricornutum. Overall, the relative abundance of psbA transcripts was greatest at 07:30 (0.133% of potential mRNA reads) (Fig. 4). Alternatively, transcripts coding for superoxide dismutase, which also catalyze the production of O 2 , were most abundant at 22:15 (Fig. 4).
Oxygen-utilization genes were diverse in type and in taxonomic distribution (Fig. 4). Transcripts for high-affinity terminal oxidases (cbb3-type cytochrome c oxidase and cytochrome d ubiquinol oxidase) as well as typical cytochrome c oxidases (COX) were transcribed at their highest levels at 12:15 (Fig. 4). As well, several monooxygenase and dioxygenase gene transcripts were detected (Fig. 4). Gene coding for the particulate methane monooxygenase, which uses O 2 in the oxidation of methane to methanol during methane oxidation (Nguyen et al., 1998;Schwarz, 2001), were most highly transcribed at 22:15 (0.718% of potential mRNA reads) (Fig. 4), and all mapped to the order Methylococcales, which includes Type I methanotrophs within the Gammaproteobacteria (Bowman, 2005). Overall, transcripts for oxygen utilization were least abundant at 07:30 (Fig. 4).

DISCUSSION
In this study, we investigated the nature and dynamics of predominant microbial taxa and metabolic processes occurring at Zodletone Spring through the analysis of mRNA transcripts expressed in sediments following sunrise (early morning), at noon, and at night. The study provides insight into a complex sulfur cycle that includes the possible involvement of Actinobacteria in sulfur transformations and identifies several interesting features not previously documented in the anoxic sediments of this high-sulfide sunlight-exposed aquatic spring, including a full methane cycle, with activity of Type I methanotrophs at night and the likelihood of O 2 production and utilization through oxygenic photosynthesis and aerobic metabolism, respectively.
Mining of metatranscriptomic datasets for gene transcripts coding for sulfur cycling processes indeed highlighted a complex sulfur cycle in Zodletone sediments. Briefly, sulfide and thiosulfate loss observed during the day (Fig. 1) could be attributed to both phototrophic and aerobic chemotrophic oxidation processes carried out by Chromatiales, Thiotrichales, and Actinobacteria (Fig. 5). Similarly, sulfidogenesis that occurs overnight ( Fig. 1) could be a function of multiple processes carried out in large part by Deltaproteobacteria, including sulfate reduction, sulfur reduction, thiosulfate disproportionation, and tetrathionate reduction (Fig. 5).
Evidence for aerobic sulfide and thiosulfate metabolism by Actinobacteria is also presented in this study (Fig. 5). Specifically, many of the transcripts encoding sulfide oxidation (flavocytochrome c sulfide dehydrogenase) and thiosulfate disproportionation (rhodanese) genes detected during the day (12:15) were mapped to the phylum Actinobacteria (Fig. 5). As taxonomic analyses of mRNA transcripts also show that Actinobacteriaaffiliated mRNA transcripts were at their highest levels during the day (Fig. 2 and Table  S3), and because many terminal cytochrome oxidase genes transcribed during the day also mapped to the phylum Actinobacteria (Fig. 4), it is also possible that these processes are supported by an aerobic mode of metabolism. Some species of Actinobacteria have been shown to be involved in aerobic oxidation of dimethylsulfide and dimethyldisulfide (Reichert et al., 1998) as well as in the oxidation of sulfur and pyrite under acidic conditions (Norris et al., 2011). Additionally, other species of Actinobacteria, Microbacterium phyllosphaerae and Leifsonia shinshuensis, were recently isolated from rhizosphere soils that are capable of the oxidation of thiosulfate, tetrathionate, trithionate and sulfur (Anandham et al., 2008). Both of these species contained a variety of genes involved in sulfur oxidation including rhodanese, thiosulfate oxidase and sulfite oxidase (Anandham et al., 2008).
These metatranscriptomic data also revealed that the production of sulfate in source sediments may not occur solely through inorganic sulfur metabolism, i.e., phototrophic sulfide/sulfur oxidation. Organo-sulfatase transcripts were much more abundant at all three time points compared to transcripts that encode for the complete oxidation of sulfur compounds to sulfate (Fig. 5). As well, organo-sulfatase transcripts were approx. twice as abundant as those that encode thiosulfate transformation (Fig. 5), another possible sulfate source. The role of sulfatases in Bacteria is typically in the cleavage of sulfate esters of organic compounds, either providing a source of sulfate or a source of organic carbon (Hanson, Best & Wong, 2004). As free sulfate is not likely a limiting nutrient in this system, the abundance of sulfatase transcripts indicates the relative importance of sulfonated organic compounds as nutrients in the spring. Together, the detection of sulfatase and sulfonation transcripts (Fig. 5) suggests that sulfate cycling includes organic sulfonated compounds within the source sediments.
Based on patterns observed from general subsystems and specific metabolic transcripts (Figs. 3A, 3B and 4) along with the shifts in taxonomic lineages observed in potential mRNA datasets (Fig. 2), we reason that the transcriptional activities of several members of the microbial community are linked to specific metabolic processes occurring during the day or night. For example, the abundance of methanogenesis functional gene transcripts (Fig. 3B) and transcripts mapping to methanogenic taxa in the early morning (Fig. 2) and the surprisingly high abundance at night of Methylococcales-affiliated transcripts and those coding for the particulate methane monooxygenase (Fig. 4), suggests that a complete methane cycle is occurring in Zodletone source sediment. Methylococcales, which has previously not been considered a dominant member in the Zodletone microbial community studies analyzing SSU rRNA genes during the daytime only (Elshahed et al., 2003;Youssef, Couger & Elshahed, 2010), may thus play a key role in the cycling of methane or other C1 compounds in Zodletone source sediments.
Despite the history of a lack of O 2 in Zodletone source sediments (Senko et al., 2004;Bühring et al., 2011), an abundance of aerobic taxa and aerobic respiration genes were detected among potential mRNA datasets (Figs. 2 and 4); these observations could be explained if oxygenic phototrophs supported the growth and metabolism of aerobic or microaerophilic taxa, a relationship that has not been previously documented under the sunlit, high-sulfide conditions that predominate at Zodletone. In this study, gene transcripts encoding oxygenic photosynthetic processes were indeed expressed by Cyanobacteria and diatoms (Fig. 4), suggesting the possibility of O 2 evolution under these conditions. While Cyanobacteria are typically considered oxygenic phototrophs, anoxygenic phototrophy has been documented among members of this phylum as an alternative to, or even replacement for, oxygenic phototrophy (Cohen, Padan & Shilo, 1975;Belkin & Padan, 1978;Jorgensen, Cohen & Revsbech, 1986;Voorhies et al., 2012). It has also been suggested that photosynthesis in Cyanobacteria will be entirely anoxygenic at sulfide concentrations >1 mM (Cohen et al., 1986). Because sulfide concentrations in Zodletone source water far exceed this threshold, ranging between 6.5-9 mM (Fig. 1), and a previous suggestion (based on the lack of detectable O 2 at the site during the daytime) that anoxygenic photosynthesis was the primary metabolic mode for Cyanobacteria in Zodletone source sediment (Bühring et al., 2011), it is surprising, here, that we report the possible production of O 2 under high-sulfide concentrations based on the transcription of cyanobacterial psbA genes, which code for the oxygen-producing subunit of Photosystem II Protein D1 (Zouni et al., 2001;Kamiya & Shen, 2003;Ferreira et al., 2004). Furthermore, it appears that this O 2 , likely produced through oxygenic photosysnthesis (Fig. 4), could be an important electron acceptor in Zodletone source sediment (Fig. 4), in spite of the highly reduced conditions present, analogous to the recently described microorganisms that can produce O 2 through nitrite reduction and then use the O 2 for methane oxidation under anoxic conditions (Ettwig et al., 2010). Because oxygen levels in the source sediment are below detectable levels, these data could suggest a yet-unidentified mechanism for oxygen transfer between oxygenic phototrophs and aerobic chemotrophs or a rapid consumption of produced O 2 , preventing its accumulation to a detectable level. While the former has not been documented in the literature, the latter has been observed in marine sediments, but with maximum rates of O 2 production exceeding those of O 2 consumption in photic zones, resulting in low steady-state O 2 concentrations (Revsbech, Madsen & Jorgensen, 1986) rather than O 2 being below detectable levels, as seen in Zodletone sediments (Senko et al., 2004;Bühring et al., 2011). However, the high sulfide concentrations at Zodletone could both greatly limit the levels of oxygenic photosynthesis to far smaller rates than observed in marine sediments or react abiotically with O 2 produced to form pentasulfide or thiosulfate, as observed in an alkaliphilic bacterial sulfoxidizing consortium (González-Sánchez & Revah, 2007).
That the O 2 available for sulfide oxidation may be derived from oxygenic photosynthesis rather than from the atmosphere or sub-oxic waters overlaying source sediments is supported by previous studies at Zodletone that find that rates of sulfide loss are most rapid in light-exposed sediments, even in aerobic incubations exposed to O 2 (Elshahed et al., 2003;Senko et al., 2004). Previous interpretations of these studies were that sulfide oxidation was driven by anoxygenic photosynthesis; the data presented here, however, also lends support to the possibility that sulfide oxidation in sunlight-exposed high sulfide environments could also be an aerobic chemolithotrophic process supported by O 2 produced from oxygenic phototrophs.

CONCLUSION
The use of metatranscriptomics for describing complex microbial ecosystems such as Zodletone spring has its limitations, including the probable occurrence of errors in functional and taxonomic assignment of functional genes. Also, the assumption involved in correlating transcript abundance with protein abundance/metabolic activity must always be taken with caution, as transcriptional and translational activity can differ in some cases (Gygi et al., 1999;Suzuki, Simon & Slabas, 2006). Despite these limitations, however, this study implicates several interesting microbial processes and interactions occurring in this sunlit, shallow, sulfur spring system, including a full sulfur and methane cycle that involve both aerobic and anaerobic bacteria in an anoxic sediment system, with O 2 that is likely produced by Cyanobacteria and diatoms by oxygenic photosynthesis. Lastly, the production of sulfate from mechanisms other than anoxygenic photosynthesis, such as aerobic sulfide and thiosulfate metabolism as well as the cleavage of sulfate esters by sulfatases, sheds insight into the microbial metabolism and lifestyles on early earth, prior to the oxygenation of the atmosphere. This work may present alternative explanations for the increase in sulfate concentrations that stimulated sulfate reduction as oxygen began to first accumulate in the atmosphere during the late Archaean and Proterozoic eons (Canfield, 1998;Reinhard et al., 2009).