# SYSTEMS BIOLOGY AND ECOLOGY OF MICROBIAL MAT COMMUNITIES

EDITED BY: Martin G. Klotz, Donald A. Bryant, Jim K. Fredrickson, William P. Inskeep and Michael Kühl PUBLISHED IN: Frontiers in Microbiology

#### *Frontiers Copyright Statement*

*© Copyright 2007-2016 Frontiers Media SA. All rights reserved. All content included on this site, such as text, graphics, logos, button icons, images, video/audio clips, downloads, data compilations and software, is the property of or is licensed to Frontiers Media SA ("Frontiers") or its licensees and/or subcontractors. The copyright in the text of individual articles is the property of their respective authors, subject to a license granted to Frontiers.*

*The compilation of articles constituting this e-book, wherever published, as well as the compilation of all other content on this site, is the exclusive property of Frontiers. For the conditions for downloading and copying of e-books from Frontiers' website, please see the Terms for Website Use. If purchasing Frontiers e-books from other websites or sources, the conditions of the website concerned apply.*

*Images and graphics not forming part of user-contributed materials may not be downloaded or copied without permission.*

*Individual articles may be downloaded and reproduced in accordance with the principles of the CC-BY licence subject to any copyright or other notices. They may not be re-sold as an e-book.*

*As author or other contributor you grant a CC-BY licence to others to reproduce your articles, including any graphics and third-party materials supplied by you, in accordance with the Conditions for Website Use and subject to any copyright notices which you include in connection with your articles and materials.*

> *All copyright, and all rights therein, are protected by national and international copyright laws.*

*The above represents a summary only. For the full conditions see the Conditions for Authors and the Conditions for Website Use.*

ISSN 1664-8714 ISBN 978-2-88919-793-4 DOI 10.3389/978-2-88919-793-4

## About Frontiers

Frontiers is more than just an open-access publisher of scholarly articles: it is a pioneering approach to the world of academia, radically improving the way scholarly research is managed. The grand vision of Frontiers is a world where all people have an equal opportunity to seek, share and generate knowledge. Frontiers provides immediate and permanent online open access to all its publications, but this alone is not enough to realize our grand goals.

## Frontiers Journal Series

The Frontiers Journal Series is a multi-tier and interdisciplinary set of open-access, online journals, promising a paradigm shift from the current review, selection and dissemination processes in academic publishing. All Frontiers journals are driven by researchers for researchers; therefore, they constitute a service to the scholarly community. At the same time, the Frontiers Journal Series operates on a revolutionary invention, the tiered publishing system, initially addressing specific communities of scholars, and gradually climbing up to broader public understanding, thus serving the interests of the lay society, too.

## Dedication to Quality

Each Frontiers article is a landmark of the highest quality, thanks to genuinely collaborative interactions between authors and review editors, who include some of the world's best academicians. Research must be certified by peers before entering a stream of knowledge that may eventually reach the public - and shape society; therefore, Frontiers only applies the most rigorous and unbiased reviews.

Frontiers revolutionizes research publishing by freely delivering the most outstanding research, evaluated with no bias from both the academic and social point of view. By applying the most advanced information technologies, Frontiers is catapulting scholarly publishing into a new generation.

## What are Frontiers Research Topics?

Frontiers Research Topics are very popular trademarks of the Frontiers Journals Series: they are collections of at least ten articles, all centered on a particular subject. With their unique mix of varied contributions from Original Research to Review Articles, Frontiers Research Topics unify the most influential researchers, the latest key findings and historical advances in a hot research area! Find out more on how to host your own Frontiers Research Topic or contribute to one as an author by contacting the Frontiers Editorial Office<EMAIL_ADDRESS>
## **SYSTEMS BIOLOGY AND ECOLOGY OF MICROBIAL MAT COMMUNITIES**

Topic Editors:

**Martin G. Klotz,** Queens College of The City University of New York (CUNY), USA **Donald A. Bryant,** The Pennsylvania State University, USA **Jim K. Fredrickson,** Pacific Northwest National Laboratory, USA **William P. Inskeep,** Montana State University, USA **Michael Kühl,** University of Copenhagen, Denmark

The primary aim of this research topic was to examine principles of systems biology in different types of stratified microbial mats present in extreme environments: high-temperature chemotrophic (A), hypersaline (B, C), and high-temperature phototrophic (D) communities. These environments exhibit contrasting gradients in key environmental variables (e.g., light, temperature, oxygen) and provide tractable systems for studying the tight metabolic coupling that occurs over relevant microbial scales.

Images adapted from:

A. Inskeep WP, Jay ZJ, Herrgard MJ, Kozubal MA, Rusch DB, Tringe SG, Macur RE, Jennings RdeM, Boyd ES, Spear JR and Roberto FF (2013) Phylogenetic and functional analysis of metagenome sequence from high-temperature archaeal habitats demonstrate linkages between metabolic potential and geochemistry. Front. Microbiol. 4:95. doi: 10.3389/fmicb.2013.00095.

B. Lindemann SR, Moran JJ, Stegen JC, Renslow RS, Hutchison JR, Cole JK, Dohnalkova AC, Tremblay J, Singh K, Malfatti SA, Chen F, Tringe SG, Beyenal H and Fredrickson JK (2013) The epsomitic phototrophic microbial mat of Hot Lake, Washington: community structural responses to seasonal cycling. Front. Microbiol. 4:323. doi: 10.3389/fmicb.2013.00323

C. Lawrence Livermore National Laboratory Genomic Sciences SFA, J. Pett-Ridge

D. Kim Y-M, Nowack S, Olsen MT, Becraft ED, Wood JM, Thiel V, Klapper I, Kühl M, Fredrickson JK, Bryant DA, Ward DM and Metz TO (2015) Diel metabolomics analysis of a hot spring chlorophototrophic microbial mat leads to new hypotheses of community member metabolisms. Front. Microbiol. 6:209. doi: 10.3389/fmicb.2015.00209

Microbial mat communities consist of dense populations of microorganisms embedded in exopolymers and/or biomineralized solid phases, and are often found in mm-cm thick assemblages, which can be stratified due to environmental gradients such as light, oxygen or sulfide. Microbial mat communities are commonly observed under extreme environmental conditions, deriving energy primarily from light and/or reduced chemicals to drive autotrophic fixation of carbon dioxide. Microbial mat ecosystems are regarded as living analogues of primordial systems on Earth, and they often form perennial structures with conspicuous stratifications of microbial populations that can be studied in situ under stable conditions for many years. Consequently, microbial mat communities are ideal natural laboratories and represent excellent model systems for studying microbial community structure and function, microbial dynamics and interactions, and discovery of new microorganisms with novel metabolic pathways potentially useful in future industrial and/or medical applications. Due to their relative simplicity and organization, microbial mat communities are often excellent testing grounds for new technologies in microbiology including micro-sensor analysis, stable isotope methodology and modern genomics. Integrative studies of microbial mat communities that combine modern biogeochemical and molecular biological methods with traditional microbiology, macro-ecological approaches, and community network modeling will provide new and detailed insights regarding the systems biology of microbial mats and the complex interplay among individual populations and their physicochemical environment. These processes ultimately control the biogeochemical cycling of energy and/or nutrients in microbial systems. Similarities in microbial community function across different types of communities from highly disparate environments may provide a deeper basis for understanding microbial community dynamics and the ecological role of specific microbial populations. Approaches and concepts developed in highly-constrained, relatively stable natural communities may also provide insights useful for studying and understanding more complex microbial communities.

**Citation:** Klotz, M. G., Bryant, D. A., Fredrickson, J. K., Inskeep, W. P., Kühl, M., eds. (2016). Systems Biology and Ecology of Microbial Mat Communities. Lausanne: Frontiers Media. doi: 10.3389/978-2-88919-793-4

# Table of Contents

Cristina Takacs-Vesbach, William P. Inskeep, Zackary J. Jay, Markus J. Herrgard, Douglas B. Rusch, Susannah G. Tringe, Mark A. Kozubal, Natsuko Hamamura, Richard E. Macur, Bruce W. Fouke, Anna-Louise Reysenbach, Timothy R. McDermott, Ryan deM. Jennings, Nicolas W. Hengartner and Gary Xie

*48 Phylogenetic and functional analysis of metagenome sequence from hightemperature archaeal habitats demonstrate linkages between metabolic potential and geochemistry*

William P. Inskeep, Zackary J. Jay, Markus J. Herrgard, Mark A. Kozubal, Douglas B. Rusch, Susannah G. Tringe, Richard E. Macur, Ryan deM. Jennings, Eric S. Boyd, John R. Spear and Francisco F. Roberto

*69 Geomicrobiology of sublacustrine thermal vents in Yellowstone Lake: geochemical controls on microbial community structure and function* William P. Inskeep, Zackary J. Jay, Richard E. Macur, Scott Clingenpeel, Aaron Tenney, David Lovalvo, Jacob P. Beam, Mark A. Kozubal, W. C. Shanks, Lisa A. Morgan, Jinjun Kan, Yuri Gorby, Shibu Yooseph and Kenneth Nealson

*85 Community structure and function of high-temperature chlorophototrophic microbial mats inhabiting diverse geothermal environments* Christian G. Klatt, William P. Inskeep, Markus J. Herrgard, Zackary J. Jay, Douglas B. Rusch, Susannah G. Tringe, M. Niki Parenteau, David M. Ward, Sarah M. Boomer, Donald A. Bryant and Scott R. Miller

*108 The epsomitic phototrophic microbial mat of Hot Lake, Washington: community structural responses to seasonal cycling* Stephen R. Lindemann, James J. Moran, James C. Stegen, Ryan S. Renslow, Janine R. Hutchison, Jessica K. Cole, Alice C. Dohnalkova, Julien Tremblay, Kanwar Singh, Stephanie A. Malfatti, Feng Chen, Susannah G. Tringe, Haluk Beyenal and James K. Fredrickson

Michael Nielsen, Niels P. Revsbech and Michael Kühl

Mohammad A. A. Al-Najjar, Alban Ramette, Michael Kühl, Waleed Hamza, Judith M. Klatt and Lubos Polerecky

Eric D. Becraft, Jason M. Wood, Douglas B. Rusch, Michael Kühl, Sheila I. Jensen, Donald A. Bryant, David W. Roberts, Frederick M. Cohan and David M. Ward

*223 The molecular dimension of microbial species: 2.* **Synechococcus** *strains representative of putative ecotypes inhabiting different depths in the Mushroom Spring microbial mat exhibit different adaptive and acclimative responses to light*

Shane Nowack, Millie T. Olsen, George A. Schaible, Eric D. Becraft, Gaozhong Shen, Isaac Klapper, Donald A. Bryant and David M. Ward

*236 The molecular dimension of microbial species: 3. Comparative genomics of*  **Synechococcus** *strains with different light responses and* **in situ** *diel transcription patterns of associated putative ecotypes in the Mushroom Spring microbial mat*

Millie T. Olsen, Shane Nowack, Jason M. Wood, Eric D. Becraft, Kurt LaButti, Anna Lipzen, Joel Martin, Wendy S. Schackwitz, Douglas B. Rusch, Frederick M. Cohan, Donald A. Bryant and David M. Ward

*249 Time dynamics of the* **Bacillus** *cereus exoproteome are shaped by cellular oxidation*

Jean-Paul Madeira, Béatrice Alpha-Bazin, Jean Armengaud and Catherine Duport

# Editorial: Systems Biology and Ecology of Microbial Mat Communities

#### Martin G. Klotz 1 \*, Donald A. Bryant 2, 3, Jim K. Fredrickson4 , William P. Inskeep5 \* and Michael Kühl 6

1 Department of Biology, Queens College of The City University of New York, New York, NY, USA, 2 Department of Biochemistry and Molecular Biology, The Pennsylvania State University, University Park, PA, USA, 3 Department of Chemistry and Biochemistry, Montana State University, Bozeman, MT, USA, 4 Pacific Northwest National Laboratory, Richland, WA, USA, 5 Department of Land Resources and Environmental Sciences, Montana State University, Bozeman, MT, USA, 6 Department of Biology, University of Copenhagen, Helsingør, Denmark

Keywords: chemotrophy, diel cycling, metabolomics, metagenomics, microbial mats, microsensors, photosynthesis, proteomics

#### **The Editorial on the Research Topic**

#### **Systems Biology and Ecology of Microbial Mat Communities**

The goals of systems biology and microbial ecology are to gain a predictive understanding of how microbial communities function in dynamic systems, inclusive of interactions among community members and spatiotemporal changes in the physicochemical environment. Microorganisms in natural systems experience cycles of environmental change over different periodicities and amplitudes, and these processes are reflected in the composition, genetic repertoire and activity of microbial community members, and community function as a whole. A primary emphasis of this research topic was to focus on reports of tractable microbial communities in extreme environments (e.g., temperature, salinity, light, pH) where a foundation of genome sequence and other molecular (-omic) and geochemical measurements provide evidence of specific functional attributes of individual community members, which are directly linked with spatiotemporal changes in key environmental variables as well as the metabolic dynamics of other community members. Thermal or saline microbial mats are often stratified with respect to key environmental variables (e.g., light, oxygen) and exhibit compositional simplicity and low heterogeneity relative to habitats such as soils and/or natural waters. Consequently, the majority of contributions to this research topic focus on either high-temperature chemotrophic microbial mats, high-temperature phototrophic mats, or hypersaline phototrophic mats in marine or epsomitic systems. The structure and function of high-temperature systems of Yellowstone National Park (YNP) was evaluated using metagenome sequence and geochemical observations across a wide range of environmental conditions, which provided a basis for understanding the distribution of thermophiles in YNP and led to the discovery of several new archaeal and bacterial lineages. Genome sequences of relevant ecotypes have provided a foundation for interrogating more detailed spatiotemporal aspects of thermophilic phototrophic communities, including microsensor analysis of their physical and chemical microenvironment, and provided a rationale for comparison to hypersaline phototrophic mats. The fixation of carbon dioxide as a primary carbon source and the production of key cofactors by autotrophs are important processes, which support diverse heterotrophs across widely different environmental circumstances. Specific metabolic linkages among community members (e.g., production of storage compounds, nitrogen fixation, fermentation, sulfate reduction, hydrogen, and vitamin production) were documented in phototrophic mats, which revealed that these

#### Edited and reviewed by:

Marc Strous, University of Calgary, Canada

#### \*Correspondence:

Martin G. Klotz<EMAIL_ADDRESS>William P. Inskeep<EMAIL_ADDRESS>
#### Specialty section:

This article was submitted to Microbial Physiology and Metabolism, a section of the journal Frontiers in Microbiology

> Received: 14 December 2015 Accepted: 22 January 2016 Published: 09 February 2016

#### Citation:

Klotz MG, Bryant DA, Fredrickson JK, Inskeep WP and Kühl M (2016) Editorial: Systems Biology and Ecology of Microbial Mat Communities. Front. Microbiol. 7:115. doi: 10.3389/fmicb.2016.00115 processes changed in unexpected ways across a diel cycle. The tight metabolic coupling among specific populations that are highly adapted to spatial and/or temporal conditions is a common theme in natural environments, and this was demonstrated in detail using different light-adapted ecotypes of cyanobacteria (Synechococcus spp.) present in alkaline siliceous geothermal mats as primary producers. Finally, controlled experiments using pure cultures and/or consortia as "systems" revealed the importance of specific nutrient requirements, and importantly, how the exoproteome of a bacterium is controlled by the level of cellular oxidation. Ultimately, a predictive understanding of the complex network of abiotic and biotic interactions that occur in natural systems will also require detailed knowledge of the regulatory and physicochemical processes that govern gene expression, posttranslational modifications, and protein activity. It is our hope that the articles included in this research topic advance a more comprehensive understanding of natural microbial communities, and demonstrate the utility of coupling molecular methods with detailed spatiotemporal measurements and dissection of microbial community function across gradients in key environmental variables such as light, temperature, pH, oxygen, or hydrogen. Integrated approaches across relevant microbial scales will lead to predictive capabilities useful for engineering microbial communities (or consortia) and for understanding how natural systems may respond to changes in key environmental variables (e.g., climate change).

## AUTHOR CONTRIBUTIONS

WI drafted the manuscript, DB, JF, MK, and MGK revised the draft and all authors agreed to the final version. The articles in the RT were edited by MGK (10), WI (3), MK (1), and DB (1).

## ACKNOWLEDGMENTS

The authors would like to thank Steve Lindemann from Pacific Northwest National Laboratory, USA for editing one of the articles in this Research Topic.

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Copyright © 2016 Klotz, Bryant, Fredrickson, Inskeep and Kühl. This is an openaccess article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

## TheYNP metagenome project: environmental parameters responsible for microbial distribution in theYellowstone geothermal ecosystem

#### **William P. Inskeep1,2\*, Zackary J. Jay 1,2, Susannah G. Tringe3\*, Markus J. Herrgård4 , Douglas B. Rusch5 and YNP Metagenome Project Steering Committee andWorking Group Members†**

1 Department of Land Resources and Environmental Sciences, Montana State University, Bozeman MT, USA

2 Thermal Biology Institute, Montana State University, Bozeman MT, USA

3 Department of Energy, Joint Genome Institute, Walnut Creek, CA, USA

4 Novo Nordisk Foundation Center for Biosustainability, Technical University of Denmark, Hørsholm Denmark

5 Center for Genomics and Bioinformatics, Indiana University, Bloomington, IN, USA

#### **Edited by:**

Martin G. Klotz, University of North Carolina at Charlotte, USA

#### **Reviewed by:**

Ludmila Chistoserdova, University of Washington, USA Gerard Muyzer, University of Amsterdam, Netherlands

#### **\*Correspondence:**

William P. Inskeep, Department of Land Resources and Environmental Sciences, Thermal Biology Institute, Montana State University, Bozeman, MT 59717, USA.

e-mail<EMAIL_ADDRESS>Susannah G. Tringe, Department of Energy, Joint Genome Institute, 2500 Mitchell Drive, Walnut Creek, CA 94598, USA.

e-mail<EMAIL_ADDRESS>
†YNP metagenome project steering committee and working group members and their respective contributions are listed in Table S1 in Supplementary Material.

The Yellowstone geothermal complex contains over 10,000 diverse geothermal features that host numerous phylogenetically deeply rooted and poorly understood archaea, bacteria, and viruses. Microbial communities in high-temperature environments are generally less diverse than soil, marine, sediment, or lake habitats and therefore offer a tremendous opportunity for studying the structure and function of different model microbial communities using environmental metagenomics. One of the broader goals of this study was to establish linkages among microbial distribution, metabolic potential, and environmental variables. Twenty geochemically distinct geothermal ecosystems representing a broad spectrum of Yellowstone hot-spring environments were used for metagenomic and geochemical analysis and included approximately equal numbers of: (1) phototrophic mats, (2) "filamentous streamer" communities, and (3) archaeal-dominated sediments. The metagenomes were analyzed using a suite of complementary and integrative bioinformatic tools, including phylogenetic and functional analysis of both individual sequence reads and assemblies of predominant phylotypes.This volume identifies major environmental determinants of a large number of thermophilic microbial lineages, many of which have not been fully described in the literature nor previously cultivated to enable functional and genomic analyses. Moreover, protein family abundance comparisons and in-depth analyses of specific genes and metabolic pathways relevant to these hot-spring environments reveal hallmark signatures of metabolic capabilities that parallel the distribution of phylotypes across specific types of geochemical environments.

**Keywords: thermophiles, geochemistry, microbial interactions, microbial mats, functional genomics**

## **INTRODUCTION**

The Yellowstone hotspot is responsible for an enormous number (>14,000) and diversity of thermal features that cover a wide range in pH (2–10), temperature (40–92˚C), and geochemical properties (Fournier, 1989; Rye and Truesdell, 2007). The waters, rocks, and mineral surfaces in these geothermal sites provide an assortment of electron donors such as hydrogen, sulfide, and ferrous iron, as well as electron acceptors (e.g., dissolved oxygen or elemental S) that are vital to the survival of thermophilic microorganisms (Brock, 1978). Remarkably, many sites exhibit relative stability or regular patterns of temporal change, despite the dynamic nature of geothermal activity (Ball et al., 1998, 2002, 2010; McCleskey et al., 2005). Microbial communities in high-temperature environments are often dominated by just several types of microorganisms (phylotypes), and in general are significantly less diverse than lowertemperature habitats. Geothermal systems thus serve as valuable models for establishing linkages between genomic potential and environmental parameters, and for understanding environmental

controls on community structure and function. Moreover, many current thermal environments inYNP harbor modern-day analogs of microbial communities potentially important in ancient Earth, perhaps most notably, cyanobacterial mats that may have analogs dating back to the beginning of oxygenic photosynthesis and the "Great Oxidation Event" (Canfield, 2005; Konhauser, 2006, 2009; Sessions et al., 2009). Consequently, the expansive field site of YNP provides a natural laboratory for studying microbial evolution and adaptation across uniquely protected and preserved thermal environments.

The evolutionary history of all three branches of Life (*Bacteria*, *Archaea*, and *Eucarya*) is intertwined with major geological events that have shaped the planet (Nealson and Ghiorse, 2001; Konhauser, 2006). Relationships between geological events and microbial function are subjects of considerable debate and interest across numerous disciplines, and provide a strong basis for interdisciplinary collaboration. The extant global distribution of microorganisms reveals an astounding diversity of functional

attributes, many of which are specialized to specific physical or chemical constraints. The importance of microbially mediated reactions in the geosciences has been discussed (Newman and Banfield, 2002; Reysenbach and Shock, 2002; Falkowski et al., 2008), and numerous metabolic pathways are directly involved in key global processes such as nutrient cycling (e.g., C, N, S, Fe), gas exchange (CO2 and CH4), and trace element cycling (e.g., As, Sb, Hg). However, the microbial complexity typical of temperate soil and aquatic habitats (e.g., Madigan et al., 2012) has contributed to the difficulty in understanding specific linkages between community members and material cycling (Tringe et al., 2005; DeLong et al., 2006; Rusch et al., 2007). Metagenomics of high-temperature communities in YNP provides an unparalled opportunity to study specific microbiological associations across widely different geochemistry and mineralogy under highly controlled, pseudo-steady-state conditions.

Metagenome sequencing of environmental DNA has provided a useful tool for studying microbial community structure and function in numerous environments (e.g., Tyson et al., 2004; Tringe et al., 2005; Baker et al., 2006; DeLong et al., 2006; Rusch et al., 2007; Dick et al., 2009). Prior metagenome sequencing in geothermal systems of YNP has been used to dissect phototrophic mat communities (Klatt et al., 2011), and establish linkages among geochemical processes and microbial populations of chemotrophic communities (Inskeep et al., 2010). These studies provide specific examples of what can be gained from detailed phylogenetic and functional analysis of metagenome sequence. Here we describe, characterize and compare metagenomes and associated metadata collected across 20 different geothermal sites (**Figure 1**), which were selected to cover some of the major types of geochemical and microbiological systems active in YNP. The specific objectives of this manuscript were to (i) discuss rationale for the study design that focused on specific linkages between environmental parameters and microbial community structure, (ii) demonstrate the primary geochemical and geophysical attributes that separate different high-temperature microbial habitats, and (iii) provide an overview of metagenome sequence content and high-level (TIGRFAM) functional analyses that serve to introduce subsequent studies focused on three main geobiological ecosystem types: (1) phototrophic mats; (2) *Aquificales*-rich "filamentousstreamer" communities; and (3) archaeal-dominated sediments. These studies are the outgrowth of one of the major collaborative activities of an NSF Research Coordination Network, established to promote and coordinate research focused on geothermal biology and geochemistry in YNP. Importantly, the research presented here represents the first comprehensive and coordinated non-PCR based survey of microorganisms distributed across a broad spectrum of high-temperature habitats in YNP.

## **RESULTS AND DISCUSSION**

#### **GEOCHEMISTRY AND ENVIRONMENTAL CONTEXT**

A geo-referenced geothermal database developed by the Yellowstone Center for Resources (U.S. National Park Service) contains information for numerous geothermal features including pH, temperature, electrical conductivity (EC), and site context (with image cataloging). A version of this database (∼8000 entries) is available

at the RCN YNP website1 , as well as other data on geothermal sites collected by either the U.S. Geological Survey (e.g., Ball et al., 1998, 2002, 2010; McCleskey et al., 2005) or individual researchers. The distribution of pH values across geothermal features in YNP is bimodal with peaks near pH 2.5 and 6.5 (**Figure 2**, based on *n* ∼ 7700). While sites chosen for this study (colored bars along *x*axis) strategically cover the bimodal pH ranges observed in YNP (**Figure 2**), the microbial communities inhabiting this wide range in pH are also strongly influenced by temperature, oxygen, sulfide, and other physical variables. EC values for a similar number of sites (*n* ∼ 6450) show that most geothermal systems in YNP exhibit EC values ranging from 0.1–3 dS/m (or mmho/cm) (**Figure 2**), which corresponds to ionic strengths (I) of ∼0.005–0.04 M. Compared to marine (*I* ∼ 0.7 M), and other saline environments, the geothermal habitats sampled herein have significantly lower levels of dissolved ions, which would not be expected to impart significant environmental selection toward halophilic (high salt) or alkaliphilic (high-pH) microorganisms.

The current study was designed to cover a broad range of geochemical and temperature conditions, and to represent several of the major types of geothermal systems in YNP (Figure S1 in Supplementary Material). The breadth of habitats sampled capitalizes on a range of several key "system-defining" variables that give rise

1http://www.rcn.montana.edu

to different niches for thermophilic microorganisms. These key environmental variables and resulting habitat types are shown as a *decision-tree* (**Figure 3**), which is a descriptive tool to visualize the study design (showing replication of similar sites when possible). The nodes in the decision-tree also serve to represent hypotheses of the primary environmental factors that are most important for defining niches occupied by specific phylotypes. The corresponding metadata (geochemistry and physical context) associated with each of these sites are provided as supplemental material (Table S2 in Supplementary Material), but key attributes (T, pH, sulfide, or elemental S, and physical context) of each site are included here to understand the predominant site groupings, as well as the targeted communities and/or populations the study aimed to elucidate (**Table 1**; **Figure 3**).

The separation of sites as a function of pH represents a key decision point in the study design. As the first node, pH allows for an important distinction between the effective temperature limit for colonization by known phototrophs in acidic (pH < 5) versus more neutral (pH 5–9) environments. It is well established that the upper temperature limit for oxygenic phototrophs is ∼70– 74˚C (Brock, 1978; Madigan et al., 2012), but cyanobacteria do not generally colonize habitats below pH 4.5–5. Consequently, at pH values < 5, a different phototrophic limit is established near 54– 56˚C, corresponding to the upper temperature limit of members of the Cyanidales (red algae), diatoms, or possible representatives from the mildly acidophilic purple-bacteria (Toplin et al., 2008; Madigan et al., 2012). Consequently, bacterial phototrophic mats are confined to the upper node of the decision-tree, and are not expected as dominant organisms in habitat types below pH 5. Diverse eukaryotic (algal) phototrophic mats exist in YNP (<55˚C), but these were not considered in the current study. The phototrophic sites (5, 6, 7, 15, 16, and 20) were chosen to target communities dominated by either oxygenic phototrophs (OP), filamentous anoxygenic phototrophs (FAPs), or anoxygenic phototrophs (AP). As noted, the sample from Obsidian Pool Prime (OPP\_17) was obtained from biofilm growth on a large glass plate and represents a different physical context relative to other phototrophic sites.

Both low- and high-pH sites are further separated based on the presence or absence of sulfide and/or elemental S (**Figure 3**). This is an important variable because high sulfide concentrations also imply low oxygen (hypoxia). When sulfide levels are above detection (or if significant elemental sulfur and/or thiosulfate are present), there is generally rapid abiotic consumption of oxygen by reduced S species (Nordstrom et al., 2005, 2009). This does not rule out the possibility that O2 influx from the atmosphere will influence the resulting microbial community, but systems with plentiful sulfide and/or elemental S will (i) contain little to no detectable dissolved oxygen, and (ii) contain species of S that can serve as potential electron donors or acceptors. Utilization of an O2 variable in the decision-tree does not offer the same utility; low O2 does not imply the presence of reduced S species. Some geothermal springs discharge non-sulfidic, hypoxic waters that contain other reduced constituents such as ferrous Fe (e.g., CP\_7, OSP\_8).

Two additional variables related to physical context are necessary to explain the study design and the corresponding habitat types included for metagenome analysis. Although these variables are included here as a simple *yes* or *no*, they represent important physical determinants of microbial community structure. Gradients in oxygen, light quantity, and quality, and sulfide are substantial as a function of mat depth and influence the distribution of microorganisms (Ward et al., 1998; Ramsing et al., 2000; van der Meer et al., 2005). For instance, sub-surface mat (SSMAT) layers are known to be less oxygenated and generally contain a greater abundance of FAPs relative to surface layers, which show a greater abundance of cyanobacteria (Ward et al., 2006; Klatt et al., 2011). Two of the phototrophic sites sampled in this study were from sub-surface mat positions at Mushroom Spring (MS\_15) and Fairy Geyser (FG\_16), where mat dissection was conducted to focus on novel FAP communities discovered in prior studies (Ward et al., 1998; Boomer et al., 2002), and because other

metagenome investigations were focused on the top layers of the MS\_15 phototrophic mats (Klatt et al., 2011; Liu et al., 2011). Other phototrophic samples (CP, WC, BLVA) were obtained from surface mats that occur within geothermal outflow channels, and thus flow is an implied physical variable that defines these habitats.

Hydrodynamic conditions (i.e., flow rate, turbulence) are crucial for defining the microenvironmental context of microorganisms in geothermal settings. All Aquificales "streamer" communities (DS\_9, MHS\_10, CS\_12, OSP\_14, OS\_11, BCH\_13) were obtained within the primary flow path of high-velocity (0.1– 0.3 ms−1 ) outflow channels, as well as site OSP\_8 (Fe-oxide mat), which is immediately down-gradient of OSP\_14. Members of the deeply rooted bacterial Aquificales are known for assuming filamentous morphology in turbulent environments where the discharge of hot reduced waters occurs in the presence of air or oxygenated waters (e.g., marine hydrothermal vents, Reysenbach et al., 2000). Flow and turbulence of discharged geothermal waters encourages the degassing of H2S(aq) as well as the ingassing of O2(g) (Inskeep et al., 2005; Nordstrom et al., 2005). Consequently, the high-velocity, stream-channel habitats are considerably different from geothermal pools, which accumulate sediments of various sulfides, S0 , alunite, kaolinite, and polymorphs of SiO2, depending on specific geochemical conditions such as pH, and levels of H2S, As, Sb, and Fe (Inskeep et al., 2009). Sulfidic and elemental S-rich sediments (CH\_1, NL\_2, MG\_3, CIS\_19, JCHS\_4, WS\_18) were sampled from low-flow habitats over a wide pH range (2.5–6.1) at temperatures considerably greater than the phototrophic limit (**Figure 3**). Sulfidic sites less than pH 5–6 can be further separated based on pH, an important variable influencing the distribution of specific crenarchaea (e.g., organisms within the order Sulfolobales increase in abundance with decreasing pH).

The oxidation of aqueous Fe(II) to solid-phase ferric oxyhydroxide is an exergonic reaction under most hydrothermal conditions (Amend and Shock, 2001; Amend et al., 2003; Inskeep et al., 2005), but of course requires dissolved O2. The in-channel flow environments promote equilibration with atmospheric conditions and encourage O2 influx. The acidic Fe-oxide mat (0.5–1 cm depth) sampled in Norris Geyser Basin (OSP\_8) occurs within a high-velocity channel where oxygenation results in O2(aq) concentrations of ∼ 40–60µM, which are 20–30% of theoretical saturation at this temperature and pressure. The concentration of Fe(II) does not appear in the decision-tree, because pH and the inference of oxygen (e.g., lack of sulfide) are already included as major determinants (concentrations of total Fe generally increase three orders of magnitude per pH unit decrease). Consequently, it is not the concentration of Fe *per se* that determines whether Fe(II) is oxidized, but whether the hydrodynamic or physical setting promotes oxygenation. For example, the more acidic CH\_1 site contains ∼5 times more soluble Fe(II) than

**Table 1 | Site names, abbreviations, and minimal set of environmental metadata necessary to separate different microbial community types (temperature, pH, dissolved sulfide or the presence of elemental S; Figure 3).**

The three major habitat types sampled include phototrophic mats (green), filamentous "streamer" communities (blue), and archaeal-dominated, elemental sulfur-rich or Fe-oxide sediments (yellow and red, respectively).

1Physical context refers to location of sample: Examples include: dissected sub-surface layer (under-mat) within a phototrophic mat; high-velocity channels containing filamentous "streamer" communities or down-gradient Fe-oxide mats; S-rich sediment in various geothermal pools

2Target phylum as discussed in submitted "RI-Files" and based on prior molecular, culture, microscopic and/or other investigations performed by individual investigators included in this project (Table S1 in Supplementary Material).

3The sample from Obsidian Pool Prime (OPP\_17) represents a different habitat type and was obtained from a large thermal pool after colonization of a glass plate. Physicochemical characteristics place this sample closest to "phototrophic mats" in the current study.

OSP\_8 (250 versus 45µM), but due to the consistent geothermal delivery of low levels of H2S(g) and the ubiquity of solid-phase elemental sulfur in CH\_1, the consumption of any O2 is likely driven by reduced phases of sulfur rather than by Fe(II). Consequently, Fe-oxides are not present in the thermal pool of CH\_1, although evidence of these phases was observed where the air-water interface meets rocks along the edge of the pool.

The study design captures a significant swath of hightemperature habitats in YNP. It is important to integrate the primary variables responsible for the separation of these habitat types together with analysis and interpretation of the metagenomes. The decision-tree (**Figure 3**) is also useful for identifying habitat types not included in the study, such as acidic geothermal habitats less than 55–60˚C, in which members of the Cyanidales (red algae) are commonly observed (Toplin et al., 2008), and in which moderately thermophilic heterotrophic bacteria are known to increase in abundance (Macur et al., 2004; Kozubal et al., 2012). The sample collected at Obsidian Pool Prime (OPP\_17) represents a different habitat context relative to other sites included in the study and was obtained from 2 months growth on a large (625 cm2 ) glass slide suspended in the water column at 56˚C (pH = 5.7). OPP is a large pool (∼0.5 ha) receiving mixed water inputs, and is the only sample in the study to be collected from the photic zone of a large water-body, albeit heavily influenced by geothermal inputs. Although the physical context is considerably different than other springs and pools included in the study, characteristics of OPP\_17 place it near the phototrophic sites in the decision-tree (**Figure 3**).

#### **SEQUENCING OVERVIEW AND MAJOR ASSEMBLIES**

One of the primary aims of the study was to utilize the distinct geochemical differences across these 20 sites as a basis for understanding the distribution and function of different thermophiles. Moreover, given the uncertainties in predicting *in situ* metabolism solely using thermodynamic favorability or inference based on distantly related cultured isolates, it is useful to determine the actual organisms and genes present in high-temperature habitats as a means of constraining the metabolic possibilities and focusing our efforts on specific electron transfer reactions or catabolic pathways. The total metagenome sequence (Sanger) per site ranged from 20 to 60 Mb (**Figure 4**), but the ratio of *assembled* sequence to total sequence varied considerably across these microbial communities. The total assembled sequence (i.e., the sum of total contig length) reflects the collapse of redundant sequence reads into contigs that correspond to the predominant populations present (**Table 1**). Consequently, for communities dominated by a single phylotype (e.g., MHS\_10, NL\_2), the total assembled sequence was only ∼2 Mbp, and there were few singlet reads remaining after assembly (nearly 90% of the sequences assembled into contigs in these cases). In contrast, more assembled sequence was obtained in more diverse communities with a greater number of dominant phylotypes (e.g., CP\_7, BLVA, OS\_11, WS\_18), but at lower contig coverage. Moreover, the significant number of singlet reads remaining after assembly in many sites suggested that the sequencing coverage was insufficient to allow genome assembly of potentially important members of these communities.

The predominant phylotypes present in these three sample groups include expected microorganisms that have been the focus of prior studies at these exact locations (**Table 1**), as well as novel microorganisms in both the *Bacteria* and *Archaea*. Members of the Chloroflexi, Cyanobacteria, Chlorobi, Firmicutes, Bacteroidetes, Acidobacteria, and Proteobacteria dominated the phototrophic mats (**Figure 5**). The only sulfidic phototrophic site (BLVA) included in the study was sampled twice (∼8 months apart), before and during a bloom of purple sulfur bacteria (Gamma-proteobacteria), and both samples contained a strong signature of green Chloroflexi. Each of the in-channel "streamer" communities were dominated by one of three major Aquificales lineages specialized to variations in pH and sulfide (**Figure 5**). However, there was considerable variation in the other organisms present across the six Aquificales communities. The low-pH

contigs > 10,000 kb. "Streamer" Communities: DS\_9 = yellow;

genomes (APIS, Badger et al., 2006).

Aquificales communities (DS\_9, OSP\_14) contained significant sub-populations of different archaea (separated by the amount of sulfide), as did sites with higher pH values (OS\_11, CS\_12, BCH\_13). In contrast, the MHS\_10 sample contained no archaeal sequence and was dominated by one genus within the Aquificales (*Sulfurihydrogenibium*-like).

The high-temperature sulfur sediments and Fe-oxide mat samples were dominated by archaeal sequence reads (i.e., >85%) and contributed an extensive diversity of archaeal protein families (TIGRFAMS) not currently represented in environmental metagenome data, with the exception of replicate sites included in a smaller study of YNP chemotrophic communities (Inskeep et al., 2010). The sites dominated by Sulfolobales (CH\_1, NL\_2) represented a low-pH extreme and these communities are compared

in more detail to other sulfur sediments at moderate, but slightly acidic pH values (pH 4–6), which are dominated by members of the orders Desulfurococcales and Thermoproteales (sites MG\_3, JCHS\_4, CIS\_19) (**Figure 5**). In contrast, the acidic Femats (OSP\_8, 14) contained Fe(II)-oxidizing Sulfolobales (*Metallosphaera* sp.) and novel archaeal populations observed only in the absence of significant levels of dissolved sulfide (**Table 1**; **Figure 5**).

#### **RELATIVE DISTRIBUTION OF PROTEIN FAMILIES**

To provide an independent assessment of the three major site groups, we analyzed broad functional differences among the 20 sites using relative TIGRFAM protein family abundances present in the metagenomes. TIGRFAM protein families represent a

relatively complete set of manually curated prokaryotic protein family models, which makes them suitable for the type of comparison we sought to perform. Abundance data, which were weighted by the scaffold assembly coverage, were obtained by counting protein copy numbers in assembled scaffolds for each site. Principal Component Analysis (PCA) of the TIGRFAM abundance data (**Figure 6**) showed that the first two principal components accounted for 80.1% of the functional variation among sites, and a two component plot (PC1 versus PC2) clearly separates the three major types of sites included in the study (phototrophic mats, "filamentous streamer communities," and archaeal-dominated S or Fe sediments). In particular, Component 2 appeared to measure the relative fraction of bacterial versus archaeal sequence in the metagenome. For example, although sites DS\_9 and OSP\_14 contain a dominant Aquificales population, these sites were more similar to the archaeal group because of the presence of subdominant (<30%) archaeal populations (albeit different archaeal phylotypes in each site). Although site OS\_11 also contained a major Aquificales population, the community contained significant sub-populations of novel bacterial groups that were also found in phototrophic mats (e.g., Bacteroidetes, Firmicutes), and thus fell closer to the phototrophic group. Site WS\_18 contained at

least three major archaeal populations and was included with the archaeal group, but this sediment sample also contained at least two major bacterial populations that plotted away from the main archaeal cluster in TIGRFAM PCA space (**Figure 6**).

Comparing groups, the Aquificales "streamer" communities and archaeal-rich habitats yielded broader functional diversity than phototrophic sites because the chemotrophic sites included both bacteria and archaea while the phototrophic sites were dominated by bacteria. The phototrophic mats contained very few sequences attributable to archaea or members of the Aquificales. Component 1 reflects the amount of Aquificales versus phototrophs (e.g., important for bacterial sites) as well as the amount of sequence data attributable to members of the Sulfolobales (especially relevant for archaeal sites), in which the content of Sulfolobales decreased from left to right (see PC1, **Figure 6**). Component 3 separates two extreme sites (CH\_1 and NL\_2) that were dominated by Sulfolobales (**Figure 5**) from the majority of other sites, although PC3 explains only 5.9% of the variation across sites.

To identify what functional categories underlie the differences among sites, we obtained average TIGRFAM (Selengut et al., 2007) abundances for each functional category and site, and clustered these values using two-way hierarchical clustering (**Figure 7**).

The site clustering tree based on TIGRFAM abundances showed remarkable similarity to the decision-tree (**Figure 3**), which indicated that the environmental variables separating the major habitat types (especially pH, temperature, and reduced sulfur versus oxygen) are key determinants of microbial community structure and function. The TIGRFAM categories that varied most significantly among sites included a broad array of proteins important in central metabolism, cell replication, motility, photosynthesis, electron transport, and other metabolic functions (**Figure 7**). The observed site separation within a TIGRFAM category could result from several factors that would not be evident without further investigation. These include the possibility that house-keeping functions simply reflect distinct phylogenetic differences observed across sites, and/or alternatively, that the range in abundance within a category may reflect the presence versus absence of specific metabolic capabilities contributed by phyla unique to a site. Consequently, the separation of sites using broad functional categories (**Figure 7**) accounts for both phylogenetic and functional differences, and provides strong support for the major site groupings emphasized in the three accompanying articles (Inskeep et al., 2013a,b; Takacs-Vesbach et al., 2013).

Specific TIGRFAM categories that varied in abundance across site groups included differences in functional categories such as nitrogen fixation and photosynthesis, expected to be highly represented in phototrophic sites (**Figure 7**). In contrast, TIGRFAMS that included cofactor biosynthesis (folic acid and lipoate), surface structures, fatty-acid biosynthesis, and chemotaxis/motility were more abundant in"streamer"communities. The archaeal sites exhibited a greater abundance of TIGRFAMs that included RNA processing, amino acid biosynthesis, nitrogen metabolism, aerobic respiration, and detoxification (**Figure 7**). In some cases, differences in TIGRFAM abundance across sites may have resulted from phylogenetic differences; for example, archaeal genes encoding a specific function may not be recognized as part of a TIGR-FAM category that has been established primarily from bacterial genomes. TIGRFAM categories that result in the greatest separation among the 20 sites included signal transduction, nitrogen fixation and regulatory functions, which represent a greater proportion of total sequences in phototrophic sites versus either Aquificales or archaeal communities (**Figure 7**; Figure S2 in Supplementary Material shows two of these TIGRFAM categories in greater detail).

More fine-grained functional differences among sites were obtained by only considering TIGRFAMs in a specific category such as "Electron Transport" (**Figure 8**). The site clustering based on the relative abundance of different electron transport domains provided a different view of the variation in functional attributes across sites. The distribution of specific respiratory complexes across sites [e.g., heme Cu oxidases (HCO), cytochrome *bd*-ubiquinol oxidases, blue Cu proteins, NiFe-hydrogenases, and nitrite/nitrate reductases] correlates more closely with geochemical parameters (e.g., oxygen, sulfur, hydrogen). Consequently, the distribution of different electron transport proteins across sites did not result in identical site clustering obtained using broader TIGRFAM categories (**Figure 7**). These observations are consistent with the fact that each group (i.e., phototrophic, archaeal, streamer communities) contained a range of sites with variable levels of sulfide or oxygen. Consequently, the factors contributing to functional variation within each of the three main site groups require additional dissection to appreciate how changes in these attributes are correlated with specific phylotypes (Inskeep et al., 2013a,b; Takacs-Vesbach et al., 2013).

The distribution and breadth of protein families observed in the current study represents a significant contribution to total protein diversity observed in all metagenomes currently found in the Integrated Microbial Genomes and Metagenomes (IMG/M) database (Markowitz et al., 2012) (as of spring 2012). This is due primarily to the fact that the current study focused on hightemperature communities rich in Aquificales, other deeply rooted novel bacteria, and numerous lineages of archaea. The functional diversity of the 20 Yellowstone geothermal samples included in this study was compared to other types of microbiomes as represented by the metagenomic sequence sets available in IMG/M, and based upon the TIGRFAM protein family abundance profiles for the assembled metagenomes. YNP samples are clearly separated from soil, water, gut, and other types of microbial communities in the public database (Figure S3 in Supplementary Material). As might be expected, the archaeal-dominated sites are most distant from the majority of public metagenomes, as few metagenome studies have targeted communities with significant archaeal populations. For example, the Crater Hills (CH\_1) and Nymph Lake (NL\_2) sites contained 1–2 predominant Sulfolobales populations and represented one extreme in the TIGRFAM PCA plots. Moreover, nearly identical samples obtained ∼2 years prior to this study (Inskeep et al., 2010) from four of the sites discussed here (CH\_1, JCHS\_4, CS\_12, MHS\_10) grouped with their expected replicate using PCA of TIGRFAM categories. Phototrophic mats and Aquificales "streamer" communities, which are dominated by eubacteria, fall closer to publicly available environmental metagenomes. The first principal component (which accounted for 24.8% of total variability in the TIGRFAM protein family abundance data across all sites in the database) almost entirely represents variation that exists within the YNP samples alone. The dramatic segregation of TIGRFAMS from YNP metagenomes versus previously sequenced mesophilic community metagenomes is indicative of the extensive contributions attributable to previously undescribed functional diversity.

## **SUMMARY**

The metagenomes and corresponding metadata discussed here and in the accompanying articles provides a significant foundation of phylogenetic and metabolic information relevant to a wide range of geothermal ecosystems in YNP. Although the total sequence data obtained on an individual site basis was not sufficient to achieve assembly of all major phylotypes present in these distinct habitats, this study makes a considerable step toward understanding how microbial community structure and metabolic potential vary across a wide range of environmental parameters. A comparative distribution of protein families identified in the metagenome sequence reflects the major differences in phylogenetic structure among three primary groups of sites studied: phototrophic mats, Aquificales "streamer"-communities, and archaeal-dominated sediments. The systematic selection of geochemically distinct sites provides an enormous opportunity to link individual phylotypes (and associated metabolic attributes) with specific physicochemical properties of the habitats they occupy.

#### **MATERIALS AND METHODS GEOCHEMICAL ANALYSIS**

Twenty geothermal sites in Yellowstone National Park (**Table 1**) were sampled and characterized in 2007–2008 (Figure S1 in Supplementary Material). Parallel samples of the bulk aqueous phase (<0.2µm) and sediment intimately associated with the microbial community were obtained simultaneously and analyzed using a

combination of field and laboratory methods. Temperature, pH, and redox-sensitive species (Fe2+/Fe3+; total dissolved sulfide; dissolved O2) were determined using field methods (with some exceptions given the number and diversity of sites studied) as described in more detail in previous reports (Langner et al., 2001; Macur et al., 2004). Total dissolved ions were determined using inductively coupled plasma (ICP) spectrometry and ion

chromatography (for all major cations, anions and trace elements). Dissolved gases (CO2, H2, CH4) were determined using closed headspace gas chromatography (Inskeep et al., 2005) of sealed serum bottles (<0.2µm) obtained in the field without air contact. A complete dataset of geochemical information corresponding to these samples is provided (Table S2 in Supplementary Material). Selected sediment and microbial mat samples were analyzed using scanning electron microscopy (Phillips Field Emission-SEM) in combination with energy-dispersive analysis of x-rays (EDAX) as well as x-ray diffraction (XRD). Thermodynamic calculations adjusted for temperature effects and performed using site-specific activities of dissolved and solid-phase constituents (Amend and Shock, 2001) showed that numerous oxidation-reduction reactions are exergonic (i.e., energy-yielding) in these types of geothermal systems and could support chemolitho- or chemoorganotrophic metabolisms (Amend et al., 2003; Inskeep et al., 2005; Shock et al., 2005). The oxidation of energy-rich, reduced constituents such as H2, CH4, H2S, S0 , and As(III) is extremely favorable when geothermal waters are exposed to atmospheric O2, or in some cases using alternate electron acceptors such as nitrate, ferric Fe, sulfate, or elemental S0 . Dissolved Fe(II) and ammonium (NH4) concentrations can be very high in certain YNP geothermal systems, and offer yet another set of exergonic reactions that serve as a potential geochemical niche for chemotrophic organisms.

#### **DNA EXTRACTION AND LIBRARY CONSTRUCTION**

A standard DNA extraction protocol was used for the majority of samples, however, several sulfur sediments were ultimately subjected to various extraction kits to generate sufficient DNA yields for library construction. Our main emphasis was to obtain representative, unbiased environmental DNA for construction of smallinsert libraries. Briefly, 3–25 g wet samples were extracted with 1 ml of Buffer A (200 mM Tris, pH 8; 50 mM EDTA; 200 mM NaCl; 2 mM sodium citrate; 10 mM CaCl2) with lysozyme (1 mg/ml final concentration) for 1.5 h at 37˚C. Proteinase K (final concentration 1 mg/ml) and SDS [final concentration 0.3% (w/v)] were then added and incubated for 0.5 h at 37˚C. This first lysate was removed and the samples were re-extracted using bead-beating protocols. The two lysates were combined and extracted with phenol-chloroform, and the resulting DNA was re-precipitated in ethanol, treated with RNAase and quantified by gel electrophoresis and staining. Our intent was to avoid biasing the samples against organisms that may be difficult to lyse with chemical methods. The extraction kits (e.g., MoBio) used on several of the sulfurcontaining sediment samples also included a modest physical lysis step. Across all 20 sites, the extracted DNA ranged in size from ∼4 to 12 kbp, and for many samples these procedures resulted in higher MW DNA ranging from 9 to 20 kbp.

#### **RANDOM SHOTGUN SEQUENCING**

Short-insert (∼3 kbp) pUC18 libraries were constructed from all sites, and preliminary sequencing of 5–10 megabases (Sanger), was performed as a quality control to determine if the initial sequencing results [via MEGAN (Huson et al., 2007) and blastx analysis] were consistent with the sample origin. All samples passed this quality checkpoint and were further sequenced to produce a total of 40–60 Mbp per site (**Figure 4**), with the exception of two sites. Both MHS\_10 and NL\_2 contained only one major population type that was covered sufficiently with ∼20 Mbp of Sanger sequence and these sites were also included for a half plate of 454 titanium pyro-sequencing (∼225 Mbp per sample), along with two additional samples (MG\_3, JCHS\_4). When combined with the total amount of Sanger sequencing performed (871 Mbp), this collectively represented nearly 2 Gbp of random shotgun sequence. 16S rRNA genes were also amplified from the DNA of each site using universal primers specific for bacteria and archaea, and one 384-well plate was sequenced from each successful library. Results from 16S rRNA gene sequencing were generally consistent with the phylogenetic signatures observed in the metagenome data, with the exception of one site (NL\_2), presumably due to PCR amplification bias. For NL\_2, the majority of 16S rRNA clones were Thermoproteales-like, but the random shotgun sequence data showed a predominant Sulfolobales-like population.

#### **SEQUENCE ASSEMBLY**

Metagenome assembly was conducted using two approaches (Celera and PGA), which resulted in reasonably similar overall assembly statistics. Due to the slightly better scaffold construction and prior history using the Celera assembler (Rusch et al., 2007; Inskeep et al., 2010), these assemblies were subjected to more detailed phylogenetic and functional analyses. However, the presence of genes or protein families described here are found in both assemblies, and the major results discussed here are not significantly different using either set of assembled data. Automated tools in IMG were used for gene identification and annotation for both sets of assembled data, and both sets of data are available on IMG/M. The following parameters were used for the Celera Assembler (Version 4.0): doOverlapTrimming = 0, doFragmentCorrection = 0, globalErrorRate = 12, utgErrorRate = 150, utgBubblePopping = 1, and useBogUnitig = 0. For PGA assemblies (Zhao et al., 2008), the following parameters were employed: OverlapLen = 30; Percent = 0.75; Clearance = 30; ClipIdn = 77; ClipQual = 10; CutoffScore = 400; EndOverhang = 800; InOverhang = 500; Min-CovRep = 50; MinLinks = 2; MinSat = 3; NumIter = 50; PenalizeN = 1; QualOverLim = 400; QualScoreCutoff = 200; Qual-SumLim = 3500; SimDiFac = 30; Verbosity = 1.

As mentioned above, four sites (NL\_2, MG\_3, JCHS\_4, MHS\_10) received both Sanger sequence and one-half plate of 454 titanium pyro-sequencing (also available on IMG/M). Due to the success of Sanger sequencing in generating significant assemblies from these four sites, the additional pyro-sequence is not discussed in great detail in the subsequent manuscripts. In all four of these sites, contig coverage increased significantly when the pyro-sequence reads were included. However, the analyses discussed here are mostly derived from the Sanger sequence data, in order to maintain a consistent analytical approach for the data from all 20 sites.

#### **PHYLOGENETIC ANALYSIS OF ENVIRONMENTAL SEQUENCE DATA**

Phylogenetic analysis of random shotgun sequence reads was performed using a number of different approaches, and these yielded convergent information regarding the predominant phylotypes present in these geothermal sites. Individual sequence reads were analyzed using G + C content (mol%) distribution coupled with blastx and MEGAN (Huson et al., 2007) analysis. Genome-level phylogenetic analysis was accomplished using fragment recruitment of environmental sequence data to reference microbial genomes (Rusch et al., 2007). At the time of writing, the database contained reference microbial genomes for ∼1500 bacteria and 100 archaea, however, only a handful of microbial genomes currently serve as appropriate references for the indigenous organisms within these communities. Assembled metagenome sequences were also analyzed using three dimensional PCA plots of nucleotide word frequencies (Teeling et al., 2004) with a simultaneous phylogenetic classification based on APIS or on a blast-based classification (Badger et al., 2006; Rusch et al., 2007). Metagenome sequence reported here can be viewed with these utilities at http://gos.jcvi.org/openAccess/scatterPlotViewer.html.

#### **HIGH-LEVEL TIGRFAM ANALYSIS OF ENVIRONMENTAL SEQUENCE DATA**

Assembled scaffolds were annotated as described in Inskeep et al. (2010) and predicted proteins from the scaffolds were assigned to TIGRFAM protein families (Selengut et al., 2007) using HMMER 3 (Eddy, 2011) with *e*-value cutoff of 1*e*−6. TIGRFAM family counts for each scaffold were multiplied by the average coverage for the scaffold and the coverage weighted counts across all scaffolds for a particular site were summed to obtain estimated total family counts for a site. Different sites were normalized so that the total count over all TIGRFAM protein families was kept constant. TIGRFAM families were categorized using a two-level functional classification described on the TIGRFAM website2 . PCA and statistical analysis of site group differences was performed using the STAMP v2.0 software (Parks and Beiko, 2010). The ANOVA test and Benjamini-Hochberg FDR correction implemented in STAMP was used to test for differences between multiple site groups. Two-way clustering was done on row-standardized (across sites) average TIGRFAM category abundance data using the Euclidean distance metric and complete-linkage hierarchical clustering. The MeV 4.8 (Saeed et al., 2003) software package was used for clustering and visualization. For comparison with other metagenomes in the IMG system, TIGRFAM profiles were generated for all the sites in this study and other metagenomic data sets downloaded using the IMG web interface.

#### **DETAILED FUNCTIONAL ANALYSIS OF ENVIRONMENTAL METAGENOMIC SEQUENCE DATA**

The assembled metagenome sequence data was also screened for specific functional genes corresponding to known or putative pathways in material and energy transfer. We were specifically interested in assessing metabolic potential for chemolithoautotrophy (CO2 fixation and electron transfer) in high-temperature geothermal systems. Query DNA sequences known to code for proteins important in the oxidation of reduced chemical constituents or the reduction of a terminal acceptor (Table S3 in Supplementary Material) were used to search the environmental sequence data. Environmental sequence fragments exhibiting sequence similarity (*e*-values < 10−10) to query sequences were then reanalyzed using blastp, and assessed individually using phylogenetic analysis of deduced protein sequences against known relatives, as well as fragment length relative to query length. False positives were minimized using this screening process. This included (i) sequences matching the correct protein family of the query sequence, but not the exact query sequence (e.g., Mo-pterin oxidoreductases versus a specific protein within this family); (ii) sequences matching a query sequence due to regions of sequence similarity, but were clearly associated with a gene or gene cluster with different function; and (iii) sequences that returned misannotated blastp relatives. It is also possible that our inventory of metabolic potential missed sequences related to a specific query gene. For example, some genes found in the metagenome data were of insufficient length relative to a specific query sequence (<40%) to make a definitive assignment. Moreover, the lower depth of coverage (<1×) of sub-dominant phylotypes precluded a complete functional analysis of these organisms.

#### **SEQUENCE AVAILABILITY**

All annotated metagenome sequence assemblies (Celera/PGA) are available through the DOE-JGI IMG/M website3 under IMG taxon OID numbers as follows (site order is identical to that presented in **Table 1**): Phototroph Sites [YNPSite06<PHONE_NUMBER>/2013515000), Site07<PHONE_NUMBER>/2014031006), Site15<PHONE_NUMBER>/2015219002), Site16<PHONE_NUMBER>/2016842 003), Site05<PHONE_NUMBER>/2013954000), Site20<PHONE_NUMBER>/2016 842008), and Site17<PHONE_NUMBER>/2016842005).]; Aquificales Sites [YNPSite09<PHONE_NUMBER>/2014031004), Site14<PHONE_NUMBER>/2013 954001),Site10<PHONE_NUMBER>/2015391001),Site12<PHONE_NUMBER>/20 14031005),Site11<PHONE_NUMBER>/2014031007),and Site13 (2022920 006/2013515002).]; Archaeal Sites [YNPSite01<PHONE_NUMBER>/2014 031002), Site02<PHONE_NUMBER>/2015219001,<PHONE_NUMBER>), Site03<PHONE_NUMBER>/2014031003,2016842001),Site19<PHONE_NUMBER>/2015 219000),Site04<PHONE_NUMBER>/2013843003),Site18<PHONE_NUMBER>/20 16842004), and Site08<PHONE_NUMBER>/2013515001).].

## **ACKNOWLEDGMENTS**

Authors appreciate support from the *National Science Foundation* Research Coordination Network Program (MCB 0342269), the DOE-Joint Genome Institute Community Sequencing Program (CSP 787081) as well as all individual author institutions and associated research support that together has made this study possible. The work conducted by the U.S. Department of Energy Joint Genome Institute is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02- 05CH11231. Authors appreciate the research permit focused on the YNP metagenome project (Permit No, YELL-5568, 2007- 2010), and managed by C. Hendrix and S. Guenther (Center for Resources, YNP).

#### **SUPPLEMENTARY MATERIAL**

The Supplementary Material for this article can be found online at http://www.frontiersin.org/Microbial\_Physiology\_and\_ Metabolism/10.3389/fmicb.2013.00067/abstract

**Table S1 | Contributions of NSF Research Coordination Network Steering Committee and Working Group members to theYellowstone Metagenome Community Sequencing Project (DOE\_JGI CSP 787081).**

**Table S2 | Geochemical parameters measured in the bulk aqueous (**<**0.2** µ**m) phase in 20 different geothermal systems sampled for metagenome analysis.**

**Table S3 | List of gene sequences and corresponding accession numbers used to query the assembled environmental sequence data for assessing**

2http://www.jcvi.org/cgi-bin/tigrfams/index.cgi

3http://img.jgi.doe.gov/m

**potential metabolic attributes associated with the predominant phylotypes found within these geothermal sites.**

**Figure S1 | Additional site photographs emphasizing landscape context of geothermal habitats and field sampling efforts (included as a separate file containing 53 annotated photographs).**

**Figure S2 |Top rankedTIGRFAM categories based on their ability to differentiate the three types of sites discussed in this study: aquificales "streamer" communities (blue), archaeal-dominated sediments (yellow),**

#### **REFERENCES**

Molecular characterization of novel red green nonsulfur bacteria from five distinct hot spring communities in Yellowstone National Park. *Appl. Environ. Microbiol.* 68, 346–355.

**and phototrophic mats (green). (A)** Signal transduction; **(B)** Regulatory functions (ANOVA test with Benjamini-Hochberg FDR multiple testing correction).

**Figure S3 | Principal components analysis of normalizedTIGRFAM protein family abundance data across all 20 sites as well as a subset of diverse metagenomic datasets from IMG.** The metagenomes are colored based on rough classification of the environments from which the sample originates: Yellowstone sites in this study = red; Yellowstone sites from other studies = orange; sediments = turquoise; soils = pink; aquatic = light-blue; fungus-garden = green; animal-gut = brown; enrichment community = light-orange; other samples = gray.

the energetics of chemolithotrophy in nonequilibrium systems: case studies of geothermal springs in Yellowstone National Park. *Geobiology* 3, 297–317.

*Front. Microbiol.* 3:109. doi:10.3389/fmicb.2012.00109

National Park," in *Geothermal Biology and Geochemistry of Yellowstone National Park*, eds W. P. Inskeep, and T. R. McDermott (Bozeman, MT: Thermal Biology Institute, Montana State University), 73–94.

*Tectonic, and Hydrothermal Processes in the YellowstoneGeoecosystem* ed. L. A. Morgan (U.S. Geol. Surv. Prof. Pap), 1717, 235–270.

tetranucleotide frequencies for the assignment of genomic fragments. *Environ. Microbiol.* 6, 938–947.

A natural view of microbial biodiversity within hot spring cyanobacterial mat communities. *Microbiol. Mol. Biol. Rev.* 62, 1353–1370.

Zhao, F., Zhao, F., Li, T., and Bryant, D. A. (2008). A new pheromone trail-based genetic algorithm for comparative genome assembly. *Nucleic Acids Res.* 36, 3455–3462.

**Conflict of Interest Statement:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

*Received: 01 December 2012; paper pending published: 06 February 2013; accepted: 09 March 2013; published online: 06 May 2013.*

*Citation: Inskeep WP, Jay ZJ, Tringe SG, Herrgård MJ, Rusch DB and YNP Metagenome Project Steering Committee and Working Group Members (2013) The YNP metagenome project: environmental parameters responsible for microbial distribution in the Yellowstone geothermal ecosystem. Front. Microbiol. 4:67. doi:10.3389/fmicb.2013.00067*

*This article was submitted to Frontiers in Microbial Physiology and Metabolism, a specialty of Frontiers in Microbiology.*

*Copyright © 2013 Inskeep, Jay, Tringe, Herrgård, Rusch and YNP Metagenome Project Steering Committee and Working Group Members. This is an openaccess article distributed under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in other forums, provided the original authors and source are credited and subject to any copyright notices concerning any third-party graphics etc.*

## Metagenome sequence analysis of filamentous microbial communities obtained from geochemically distinct geothermal channels reveals specialization of three Aquificales lineages

**Cristina Takacs-Vesbach1\*,William P. Inskeep2\*, Zackary J. Jay 2 , Markus J. Herrgard3 , Douglas B. Rusch4 , Susannah G. Tringe5 , Mark A. Kozubal 2 , Natsuko Hamamura6 , Richard E. Macur 2 , BruceW. Fouke7 , Anna-Louise Reysenbach8 ,Timothy R. McDermott 2 , Ryan deM. Jennings 2 , NicolasW. Hengartner 9 and Gary Xie10**

1 Department of Biology, University of New Mexico, Albuquerque, NM, USA

5 Department of Energy-Joint Genome Institute, Walnut Creek, CA, USA

6 Center for Marine Environmental Studies, Ehime University, Matsuyama, Ehime, Japan

7 Roy J. Carver Biotechnology Center, University of Illinois, Urbana, IL, USA

8 Department of Biology, Portland State University, Portland, OR, USA

9 Computer, Computational, and Statistical Sciences Division, Los Alamos National Laboratory, Los Alamos, NM, USA

10 Bioscience Division, Los Alamos National Laboratory, Los Alamos, NM, USA

#### **Edited by:**

Martin G. Klotz, University of North Carolina at Charlotte, USA

#### **Reviewed by:**

Olivia Mason, Florida State University, USA

Martin Keller, Oak Ridge National Laboratory, USA

#### **\*Correspondence:**

Cristina Takacs-Vesbach, Department of Biology, University of New Mexico, MSC03 2020, Albuquerque, NM 87131, USA. e-mail<EMAIL_ADDRESS>William P. Inskeep, Department of Land Resources and Environmental Sciences, Thermal Biology Institute, Montana State University, Bozeman, MT 59717, USA.

e-mail<EMAIL_ADDRESS>
The Aquificales are thermophilic microorganisms that inhabit hydrothermal systems worldwide and are considered one of the earliest lineages of the domain Bacteria. We analyzed metagenome sequence obtained from six thermal "filamentous streamer" communities (∼40 Mbp per site), which targeted three different groups of Aquificales found in Yellowstone National Park (YNP). Unassembled metagenome sequence and PCR-amplified 16S rRNA gene libraries revealed that acidic, sulfidic sites were dominated by Hydrogenobaculum (Aquificaceae) populations, whereas the circum-neutral pH (6.5–7.8) sites containing dissolved sulfide were dominated by Sulfurihydrogenibium spp. (Hydrogenothermaceae). Thermocrinis (Aquificaceae) populations were found primarily in the circum-neutral sites with undetectable sulfide, and to a lesser extent in one sulfidic system at pH 8. Phylogenetic analysis of assembled sequence containing 16S rRNA genes as well as conserved protein-encoding genes revealed that the composition and function of these communities varied across geochemical conditions. Each Aquificales lineage contained genes for CO2 fixation by the reverse-TCA cycle, but only the Sulfurihydrogenibium populations perform citrate cleavage using ATP citrate lyase (Acl). The Aquificaceae populations use an alternative pathway catalyzed by two separate enzymes, citryl-CoA synthetase (Ccs), and citryl-CoA lyase (Ccl). All three Aquificales lineages contained evidence of aerobic respiration, albeit due to completely different types of heme Cu oxidases (subunit I) involved in oxygen reduction. The distribution of Aquificales populations and differences among functional genes involved in energy generation and electron transport is consistent with the hypothesis that geochemical parameters (e.g., pH, sulfide, H2, O2) have resulted in niche specialization among members of the Aquificales.

**Keywords: thermophiles, functional genomics, phylogeny, autotrophic processes, sulfide oxidation**

#### **INTRODUCTION**

The order Aquificales represents a group of thermophilic microorganisms that inhabit marine and terrestrial hydrothermal systems worldwide (Ferrera et al., 2007). This lineage is of significant interest because its members are believed to comprise the deepest lineage of the domain *Bacteria* (Coenye and Vandamme, 2004; Barion et al., 2007), although alternative evolutionary histories have been suggested (Griffiths and Gupta, 2004; Boussau et al., 2008; Zhaxybayeva et al., 2009). The Aquificales include two predominant families, the Hydrogenothermaceae and Aquificaceae, and both are well-represented in different geothermal features of Yellowstone National Park (YNP) (Reysenbach et al., 2005). Members of the Hydrogenothermaceae include the genus *Sulfurihydrogenibium*, which inhabit circum-neutral sulfidic springs in YNP (Hugenholtz et al., 1998; Reysenbach et al., 2000b, 2005). The Aquificaceae comprise two divergent groups:

*Hydrogenobaculum* spp. are found predominantly in low-pH systems (pH < 4) (Jackson et al., 2001; Langner et al., 2001; Macur et al., 2004; D'Imperio et al., 2007, 2008; Hamamura et al., 2009), while *Thermocrinis-*like organisms generally exhibit higher-pH ranges (pH 6–9), overlapping with the near circum-neutral optimum for *Sulfurihydrogenibium* (Reysenbach et al., 1994; Blank et al., 2002).

Given our limited understanding of their phylogenetic history and metabolic potential, numerous questions remain regarding the ecology and evolution of the Aquificales. For example, the Aquificales are important primary producers in hydrothermal systems (Harmsen et al., 1997; Yamamoto et al., 1998; Reysenbach et al., 1999, 2000a; Blank et al., 2002; Eder and Huber, 2002; Inagaki et al., 2003; Spear et al., 2005) and several members of this lineage have been shown to use the reductive tricarboxylic acid (r-TCA) cycle for the fixation of carbon dioxide (Shiba et al., 1985; Beh et al., 1993; Ferrera et al., 2007). However, several members of the order are also capable of heterotrophy (Huber et al., 1998; Nakagawa et al., 2005; Caldwell et al., 2010), which makes their ecological role as possible primary producers unclear. Moreover, the diversity and variation of r-TCA-specific enzymes and substrates used by cultured members has led to speculation regarding the origin and the distribution of this pathway within the phylum (Hugler et al., 2007).

The Aquificales were originally named for the ability of the type strain, *Aquifex pyrophilus,* to oxidize molecular hydrogen to water (Huber et al., 1992). Consequently, hydrogen oxidation was generalized to the entire order based on the phenotype of a few cultivated members (Reysenbach and Cady, 2001; Donahoe-Christiansen et al., 2004; Huber and Eder, 2006; D'Imperio et al., 2008). Additionally, the abundance of molecular hydrogen in hydrothermal systems was taken as evidence for the predominance of this metabolism in thermal features from YNP (Spear et al., 2005). However, the recent isolation of several new terrestrial species indicates some organisms within this group are unable to grow on hydrogen and do not contain Group I Ni-Fe hydrogenases (Reysenbach et al., 2009). Diverse metabolisms have been detected among cultured members of the Aquificales including the use of H2, elemental sulfur, and thiosulfate as energy sources, and although generally aerobic, some members are microaerophilic and/or utilize nitrate as an electron acceptor (D'Imperio et al., 2008; Reysenbach et al., 2009). Thus, while the metabolic potential within the Aquificales is largely known from a few well-studied isolates, there is not an equal distribution of representatives from each of the families in culture, and their physiological capabilities remain unknown, especially under the environmental conditions that are characteristic of their natural habitats.

Much of what is known about the distribution of the Aquificales is based on molecular diversity studies across different habitat types (Reysenbach et al., 1994, 2000b; Hugenholtz et al., 1998; Stohr et al., 2001; Van Dover et al., 2001; Inagaki et al., 2003), although their ecology has been inferred largely from the geochemical conditions they inhabit (Fouke et al., 2000; Spear et al., 2005; Hall et al., 2008; Hamamura et al., 2009) and the physiology of the few cultivated members (Jahnke et al., 2001; Takai et al., 2002; Reysenbach et al., 2009). Based on molecular diversity surveys of 16S rRNA and metabolic genes, different Aquificales lineages do not generally share the same habitats in YNP (Reysenbach et al., 2005), although exceptions have been noted where *Thermocrinis-* and *Sulfurihydrogenibium*-like organisms have been

**FIGURE 1 | Site photographs of Aquificales "streamer" communities sampled fromYellowstone National Park**. The sites represent diverse geochemical environments as noted from the different mineralogy apparent visually [DS\_9 (elemental sulfur); OSP\_14 (Fe-oxides); MHS\_10 (calcium carbonate); CS\_12 (pyrite, amorphous Fe-sulfides); OS\_11 (none); BCH\_13 (none, Fe(III)-staining on silica)].

found together in sulfidic circum-neutral pH channels (Hall et al., 2008; Hamamura et al., 2009). Therefore, the distinct separation of different Aquificales lineages as a function of geochemical conditions provides a unique opportunity to study the evolutionary and ecological history of diverse members of this group *in situ*. Here we provide a phylogenetic and functional analysis of metagenome sequence obtained from six Aquificales "streamer" communities, representing two replicate communities of each of the three major Aquificales lineages found in high-velocity outflow channels of YNP. The distribution and metabolic potential of Aquificales in these distinct habitat types was correlated with the geochemical attributes and conditions measured in these same locations.

#### **RESULTS**

#### **ENVIRONMENTAL AND GEOCHEMICAL CONTEXT**

Each high-temperature microbial community (including associated biomineralized solid-phases) was sampled from the primary flow-path (e.g., Veysey et al., 2008) of geothermal outflow channels with high velocities ranging from 0.2 to 0.5 ms−1 . The sites were chosen to obtain a range in pH (3–8) and other geochemical attributes such as dissolved oxygen, dissolved sulfide, and/or predominant solid-phases associated with each microbial community (**Figure 1**; **Table 1**). Moreover, the physical and geochemical characteristics of these sites provide a representative subset of major Aquificales habitat types common in YNP (Reysenbach et al., 2005). Filamentous Aquificales"streamer"communities are generally found along the primary flow paths of thermal channels and colonize hydrodynamic regimes of shallow (<2 cm deep), highvelocity, turbulent (Reynold's Number > 100,000; Fouke, 2011) spring water. These zones exhibit rapid outgassing of dissolved gases such as CO2 and H2S, and in-gassing of oxygen (Inskeep et al., 2005; Kandianis et al., 2008; Fouke, 2011).

The acidic (pH 3–3.5) sites in Norris Geyser Basin (NGB) represent two common Aquificales habitats in YNP. The relatively high concentration of H2S(aq) in *Dragon Spring* (DS\_9) (80– 100µM) results in the deposition of copious amounts of elemental sulfur due to both biotic and abiotic oxidation (D'Imperio et al., 2008), whereas lower concentrations of H2S(aq) found in the *One Hundred Spring Plain* (OSP\_14) site (<10µM) are not sufficient to form elemental sulfur. In contrast, the "streamer" communities at OSP\_14 form at a transition from reduced to oxygenated source waters and subsequent deposition of Fe-oxides (**Figure 1**). Scanning electron microscopy, elemental analysis, and electron diffraction confirm that elemental sulfur is the dominant solid phase in DS\_9, and amorphous Fe(III)-oxide is the predominant phase associated with the community at OSP\_14 (Langner et al., 2001; Inskeep et al., 2005) (**Figure 1**), although the gray discoloration of these Fe-oxides likely results from interaction of low concentrations of H2S with ferric oxides.

The outflow channels at *Mammoth Hot Springs* (MHS\_10) and *Calcite Springs* (CS\_12) both contain high concentrations of dissolved sulfide (>100µM), but other geochemical differences result in different solid-phases biomineralized within the"streamer"fabric (**Figure 1**). The hydrothermal fluids discharging at MHS intersect the Madison Limestone formation (Fournier, 1989; Fouke, 2011) and as a result, contain high concentrations of dissolved inorganic carbon (DIC), Ca, and to a lesser extent Mg. Upon discharge, CaCO3 is biomineralized as aragonite needles within and along microbial filaments forming the"streamer"fabric (**Figure 1**) (Fouke et al., 2003; Kandianis et al., 2008; Fouke, 2011). In contrast, the black, filamentous structures at CS\_12 are comprised primarily of pyrite formed along and within intertwined microbial filaments. The higher-pH, dissolved Fe, and total dissolved sulfide at *Calcite Spring* all favor the precipitation of iron sulfides relative to MHS\_10 where no pyrite is observed. Rhombohedral crystals of elemental sulfur were found within the "streamers" of both MHS\_10 and CS\_12; however, the precipitation of elemental sulfur due to reaction of dissolved sulfide with oxygen is favored at low-pH (Xu et al., 1998; Nordstrom et al., 2005) when H2S(aq)/HS− 1 (the pKa of H2S at 70˚C is ∼6.8; Amend and Shock, 2001).

Higher-pH (pH 7.5–8) "streamer" communities were also sampled from the outflow channels of low sulfide, alkaline siliceous

1 I, ionic strength calculated from aqueous geochemical modeling at sample temperature; DIC, dissolved inorganic C; DS, dissolved sulfide; S2O3, thiosulfate; DO, dissolved oxygen; CH4 and H2 values are for aqueous species.

2Predominant solid phases determined using scanning electron microscopy (FE-SEM) coupled with energy dispersive analysis of x-rays (EDAX) and x-ray diffraction (XRD).

springs (*Octopus* and *Bechler Springs* (OS\_11 and BCH\_13). Under these conditions, filamentous communities within the flow channels at temperatures ranging from 78 to 82˚C appear light-pink to gelatinous,and exhibit no obvious deposition of reaction products that may have formed from the oxidation of reduced species (e.g., no Fe-oxides, elemental sulfur or pyrite solid-phases are formed). This is consistent with the geochemical characteristics of these systems indicating low concentrations of reduced species such as sulfide, hydrogen and Fe(II), and mainly comprised of dissolved Na, Cl, and Si (**Table 1**). The predominant solid-phase(s) deposited along the evaporative margins of these channels are various forms of silica (although Fe-oxide staining is evident in BCH\_13).

#### **TAXONOMIC CLASSIFICATION OF METAGENOME SEQUENCE**

We first report on analysis of individual metagenome sequences (average length ∼800 bp), then summarize results from metagenome assemblies. The G + C content (%) was determined for individual Sanger sequences (∼800 bp per read) from each of the six sites and classified taxonomically using "blastx" against the NCBI database (**Figure 2**). In addition, individual sequences were compared to available reference genomes using fragmentrecruitment analysis (nucleotide) against more than 1400 reference genomes (tools developed during the Global Ocean Survey; Rusch et al., 2007) (**Figure A1** in Appendix).

Classification of unassembled DNA sequence reads (**Figure 2**) showed that the acidic sites (DS\_9 and OSP\_14) were both dominated by *Hydrogenobaculum* populations with an average G + C content (%) of 34.7%. However, these two samples were considerably different with respect to the archaeal populations present. The elemental sulfur site (DS\_9) contained novel populations within the Euryarchaeota and Thaumarchaeota, whereas the Fe-rich streamers (OSP\_14) contained a *Metallosphaera yellowstonensis*-like population known to be an important Fe(II)-oxidizer in these habitats (Kozubal et al., 2011). Two circum-neutral sites (pH 6.5 and 7.8), which both contain high concentrations of dissolved sulfide, provided an excellent comparison to the lower pH (3–3.5) sulfidic habitats. A major shift from *Hydrogenobaculum*-like organisms at pH 3–3.4 (DS\_9, OSP\_14) to *Sulfurihydrogenibium*-like organisms at pH 6.5 and 7.8 (MHS\_10 and CS\_12) was observed across these environments. Nearly 90% of the unassembled Sanger sequence reads from MHS\_10 were related to *Sulfurihydrogenibium* spp. (**Figure 2**). The G + C content distribution plots of sequences from MHS\_10 and CS\_12 show the importance of *Sulfurihydrogenibium*-like organisms with an average G + C value of 32.4%. A second major Aquificales population is evident in CS\_12 with an average G + C content of 46.5%, and is closely related to the major Aquificales organisms present in OS\_11 and BCH\_13 (**Figure 2**). The relative abundance of Aquificales lineages observed using all random sequences was generally very similar to the distribution of 16S rRNA gene sequences PCR-amplified using universal bacterial primers (**Figure A4** in Appendix).

The outflow channels of *Octopus* (OS\_11) and *Bechler* (BCH\_13) *Springs* exhibit similar pH values to CS\_12 (pH 7.5–8),

**National Park (YNP)**. Phylogenetic classification of each sequence read (∼800 bp) was performed using MEGAN (based on "blastx"), and provides one method for visualizing the predominant phylotypes that contribute to

green = Sulfurihydrogenibium; violet =Thermus; blue =Aquificaceae; light-blue = Hydrogenivirga; light-orange = Domain Bacteria; light-green =Thermoproteaceae).

but do not contain significant concentrations of dissolved sulfide (below detection). The predominant Aquificales sequences in OS\_11 and BCH\_13 have nearly identical G + C content distributions (average G + C = 46.7%) and taxonomic assignment (*Thermocrinis*-like), and were also similar to the sub-dominant Aquificales population in CS\_12 (**Figure 2**). Past work on the distribution of 16S rRNA genes in OS\_11, BCH\_13, and similar channels have also shown that these Aquificales are related to *Thermocrinis* spp. (Hugenholtz et al., 1998; Reysenbach et al., 2005;Hall et al., 2008). The only genome sequence available within this genus was *T. albus*, and this genome did not serve as an adequate reference for metagenome sequence from these sites (e.g., see fragment recruitment, **Figure A1** in Appendix). Consequently, the YNP *Thermocrinis*-like populations are not well-represented by currently available reference genomes, and phylogenetic assignment of these population(s) in the G + C content frequency plot (**Figure 2**) was made at the family level (Aquificaceae). In contrast, the genome sequence of *Hydrogenobaculum* sp. Y04AAS1 is a reasonably good reference for populations present in DS\_9 and OSP\_14, while the *Sulfurihydrogenibium* spp. genomes (strain Y03AOP1 or *S. yellowstonense;* Reysenbach et al., 2009) serve as good references for metagenome sequence from MHS\_10 and CS\_12 (**Figure A1** in Appendix). The recently released *Hydrogenobaculum* spp. genomes isolated from *Dragon Spring* (Romano et al., 2013) were not compared to DS\_9 assemblies here, but these reference genomes are likely superior (relative to strain Y04AAS1) for direct comparison to *Hydrogenobaculum*-like sequence from DS\_9.

#### **PHYLOGENETIC ANALYSIS OF METAGENOME SEQUENCE ASSEMBLIES**

The *de novo* assemblies from each Aquificales community resulted in significant consensus sequence of indigenous populations present in these environments, regardless of whether the Celera (Rusch et al., 2007) or the PGA assembler (Zhao et al., 2008) was used. A similar number and distribution of contigs was generated from the two assemblies. In several cases, the amount of consensus sequence corresponding to the major Aquificales populations approached the expected size of reference Aquificales genomes (Reysenbach et al., 2009). Ordination of nucleotide word frequency analysis (NWF-PCA, or k-mer analysis) (Teeling et al., 2004) has been shown to be an excellent criterion for determining whether contigs generated during assembly may belong to the same species, and was used to evaluate all Aquificales contigs greater than 2000 bp (**Figure 3**). Similar Aquificales populations were observed in each of the replicate field sites: DS\_9 (yellow) and OSP\_14 (red); MHS\_10 (green) and CS\_12 (violet), and OS\_11 (dark-blue) and BCH\_13 (light-blue), respectively. Moreover, clear separation in Aquificales sequence was observed across sites with major geochemical differences (**Figure 3A**).

Phylogenetic analysis of these sequence clusters using the Automated Phylogenetic Identification System (APIS; Badger et al., 2006) in identical PCA orientation shows that the majority of Aquificales contigs from sites DS\_9 and OSP\_14 were most closely related to *Hydrogenobaculum* sp. Y04AAS1 (**Figure 3B**). Sequence clusters from MHS\_10 consistently showed highest similarity to either the *Sulfurihydrogenibium* sp. Y03AOP1 or *S. yellowstonense* reference genomes. CS\_12 contained Aquificales contigs most

similar to the *Sulfurihydrogenibium* sp. genomes, as well as a separate population showing greater identity to the *Aquifex aeolicus* VF5 genome (i.e.,Aquificaceae). In the absence of a reference strain with higher sequence identity, the *Aquifex aeolicus* VF5 genome served to represent this *Thermocrinis*-like population (see also **Figure 2**).

Contigs from the higher-pH (pH ∼8), non-sulfidic flow channels (OS\_11 and BCH\_13) cluster together and were also classified at the family level, Aquificaceae (**Figure 3**). As discussed above, the Aquificales populations within OS\_11 and BCH\_13 are actually related to *Thermocrinis* spp., as determined from phylogenetic analysis of 16S rRNA and other genes present within each of the metagenome assemblies (**Figure 4**). However, there are currently no adequate reference strains for comparison to these YNP populations. For example, the recently released *T. albus* genome is not a good reference for the YNP *Thermocrinis* populations from OS\_11 or BCH\_13, based on poor nucleotide identity and lack of similar gene content (**Figure A2** in Appendix).

The composition of each community was also evaluated using phylogenetic analysis of 16S rRNA genes detected within the assembled sequence. Of 407 16S rRNA genes identified by the integrated metagenome data management and comparative analysis system (IMG/M, Markowitz et al., 2012) from all six sites, 64 were of sufficient length (>1100 bp) for phylogenetic analysis. Aquificales 16S rRNA gene sequences dominated the dataset (**Figure 4**) and grouped with *Hydrogenobaculum* (DS\_9 and OSP\_14), *Sulfurihydrogenibium* (MHS\_10 and CS\_12), or *Thermocrinis* spp. (OS\_11, BCH\_13, CS\_12). Unclassified 16S rRNA gene sequences grouped with novel DNA sequences reported in previous studies (Reysenbach et al., 1994; Blank et al., 2002; Hall et al., 2008) and represent members of potentially deeply rooted candidate phyla (Hall et al., 2008). Finally, PCR-amplified bacterial and archaeal 16S rRNA gene clone libraries largely agreed with the diversity detected within the assembled metagenome sequence (**Figures A3**

and **A4** in Appendix). Although archaeal primers failed to amplify from several samples, metagenome sequence contained significant levels of different archaea (**Figure 2**).

To broaden our phylogenetic analysis, the community composition was also investigated by assigning phylotypes of 31 conserved protein-encoding marker genes using AMPHORA (Wu and Eisen, 2008). Of the 159,636 predicted genes among the six metagenomes, 994 were protein markers detected by AMPHORA and of sufficient length for phylogenetic analysis. At least 29 of the 31 housekeeping genes included in the AMPHORA database were detected in each of the six sites. The majority of marker genes for each site were assigned to members of the Aquificales, which reflects the greater amount of assembled Aquificales-like sequence across these sites, relative to other organisms (**Figure A5** in Appendix). Significant amounts of assembled sequence corresponding to members of the Firmicutes, Alphaproteobacteria, Thermotogae, Bacteroidetes, Spirochetes, Proteobacteria, Actinobacteria, Deinococcus-Thermus, Acidobacteria, Gammaproteobacteria, Planctomycetes, Chlamydiae, Betaproteobacteria, Fusobacteria, and Epsilonproteobacteria, as well as unclassified *Bacteria* and *Archaea* was predicted using AMPHORA (**Figure A5** in Appendix). A minor number of sequence matches to Chloroflexi, Chlorobi, and/or Cyanobacteria may reflect exogenous inputs from lower temperatures. The version of AMPHORA used here suggests greater diversity than indicated by other more direct phylogenetic analyses of either individual sequences, assembled contigs, or PCR-amplified 16S rRNA genes (e.g., **Figures 2** and **3**; **Figures A2**–**A4** in Appendix), and is due in part to a limited number of appropriate reference genomes included in the AMPHORA database, especially from thermophiles and other deeply rooted, uncharacterized lineages. This prohibits accurate phylogenetic placement of a subset of these sequences. Consequently, the actual abundance of each population in these samples is better understood in terms of total sequence reads (i.e., **Figure 2**).

Population-level community richness can be inferred from the average abundance of single-copy genes detected within a specific population, assuming sufficient coverage (Table S1 in Supplementary Material). Given that a majority of the AMPHORA protein markers were detected for the Aquificales populations present in each site, the abundance of single-copy genes should provide a reasonable estimate of population-level heterogeneity. AMPHORA markers suggested that MHS\_10 is dominated by a single *Sulfurihydrogenibium* sp., while CS\_12 averaged five distinct *Sulfurihydrogenibium*-related copies of each single-copy gene (Table S1 in Supplementary Material). The acidic communities (DS\_9 and OSP\_14) may be comprised of three distinct *Hydrogenobaculum*-like populations, whereas OS\_11 and BCH\_13 contained at least eight or nine novel *Thermocrinis*-related populations. The number of Aquificales 16S rRNA genes identified in the metagenome assemblies (**Figure 4**) did not always agree with the number of Aquificales-related single-copy genes detected using AMPHORA. For example, 13 *Sulfurihydrogenibium*-like 16S rRNA genes were detected in MHS, yet the average abundance of single-copy genes suggest that this site is dominated by only one *Sulfurihydrogenibium*-like population. This discrepancy may be attributed to the fact that 16S rRNA genes do not often assemble with the same consistency as other housekeeping genes (Rusch et al., 2007), and/or to the possible presence of multiple rRNA operon copies (sequenced *Sulfurihydrogenibium* genomes have two to four annotated 16S rRNA gene copies). Although the sequencing depth obtained in the current study is not sufficient to accurately assess these communities at the ecotype level (Ward et al., 2006), it is clear that sequence variants of highly related Aquificales populations are evident within each site.

#### **FUNCTIONAL ANALYSIS OF "STREAMER" COMMUNITIES IN YNP**

The functional gene content among the Aquificales "streamer" communities was compared using multivariate statistical analysis of protein family (TIGRFAM) abundances based on all predicted proteins from assembled metagenome sequence. To study different driving forces behind functional diversity, the analysis was performed using both the full set of all TIGRFAM protein families and a set of TIGRFAMS associated with electron transport (ET) functions. Both analyses reveal contributions from the predominant Aquificales lineage(s) across sites as well as the variable co-community members observed within each individual site (e.g., **Figure 2**). Site comparisons were made using PCA of all TIGRFAMS, which show strong support for the similarities between site pairs in acidic (DS\_9, OSP\_14) and circum-neutral sulfidic systems (MHS\_10, CS\_12) (**Figure 5A**). The first principal component explained over 57% of the total variation in protein family abundance among the different "streamer" communities and separated the two predominant families of Aquificales represented across these sites (Aquificaceae include both the *Hydrogenobaculum* sp. that dominate sites DS\_9 and OSP\_14, and the *Thermocrinis*-like populations in OS\_11 and BCH\_13). Also, these four sites all contained archaeal co-community members (e.g., ∼30% of sequence reads in OSP\_14 and DS\_9, ∼20% in OS\_11, and ∼8% in BCH\_13). Conversely, sites MHS\_10 and CS\_12 were dominated by *Sulfurihydrogenibium* populations in

the family Hydrogenothermaceae. Despite the similar *Thermocrinis*-like populations present in OS\_11 and BCH\_13, these communities exhibited considerable functional differences, which ultimately tracked with differences in co-community members. Several of the phylogenetically distinct bacterial populations observed in OS\_11 and BCH\_13 (see **Figure 2**), which are not well characterized, contribute to unexpected functional differences between these communities. Factor 2, which accounted for ∼22% of variation in relative TIGRFAM abundance across these six sites, separated the low-pH *Hydrogenobaculum* sites from higher-pH sites.

Principal components analysis (PCA) using a subset of "ET" TIGRFAMS shows a slightly different pattern in site separation (**Figure 5B**) that is more sensitive to specific electron donors and acceptors (i.e., respiratory pathways), and more consistent with geochemical differences across sites (**Table 1**). Factor 1, which represented 63% of variation across sites, separated sites based on pH (low-pH sites DS\_9 and OSP\_14 versus other four sites), while Factor 2 correlated strongly with the presence or absence of sulfide, and separated the more oxic sites OS\_11 and BCH\_13 from sulfidic sites. The third principal component (∼6% of functional variation) emphasized differences in community composition between the two oxic sites (OS\_11 and BCH\_13), and did not contribute to overall functional differences in the more sulfidic systems (DS\_9, OSP\_14, MHS\_10, CS\_12).

Hierarchical cluster analysis of the "ET" TIGRFAMs (**Figure 6**) show expected site pairing based on geochemical attributes identified in the study design (Inskeep et al., 2013), and can be used to identify specific protein families that explain the overall variation seen in PCA (**Figure 5B**). These TIGRFAMs included ET and terminal oxidase proteins specific to different respiratory pathways dependent on either sulfur, arsenite, hydrogen, and oxygen as well as other cytochromes, ferrodoxins, and flavoproteins (**Figure 6**; Table S2 in Supplementary Material). For example, while MHS\_10 and CS\_12 sites were both dominated by *Sulfurihydrogenibium* sp., the secondary *Thermocrinis* (family Aquificaceae*)* and *Thermus* populations in CS\_12 contributed a significant number of additional cytochrome families that are not present in *Sulfurihydrogenibium* sp. Similarly, the higher complexity community at OS\_11 contained several additional gene families compared to the lower complexity *Thermocrinis*-dominated community at BCH\_13.

Hierarchical cluster analysis using broader TIGRFAM profiles (**Figure A6** in Appendix) supports the relative similarity of acidic sites (DS\_9 and OSP\_14) and circum-neutral sulfidic sites (MHS\_10 and CS\_12); however, the higher-pH non-sulfidic (i.e., oxic) communities (OS\_11, BCH\_13) did not form a separate group. Despite similar *Thermocrinis* populations in these two communities (e.g., see **Figure 2**), the overall functional profiles are quite different and this is consistent with the diverse and novel bacterial phylotypes observed in these sites (especially OS\_11). Detailed analysis of specific functional categories across sites reveal differences in the relative abundance of genes coding for numerous cellular functions including nucleotide and DNA metabolism, regulatory functions, energy metabolism, central C metabolism, mobile elements, transcription, cofactors, and transporters (**Figure A6** in Appendix).

The TIGRFAMs with the greatest statistical site separation (lowest *p*-value) include the pyridine nucleotide biosynthesis and PTS signal transduction categories (see **Figure A7** in Appendix and Table S3 in Supplementary Material) for a complete list of TIGRFAM *p*-values). A coarse view of the importance of specific geochemical variables such as pH can be examined by looking at specific TIGRFAMs providing the highest statistical separation of sites based on selected variables. For example, the most pH-dependent TIGRFAMs include categories such as glutathione disulfide reductases, thioredoxin-disulfide reductases, formate dehydrogenases, and NADH-plastoquinone oxidoreductases (**Figure A8** in Appendix).

#### **FUNCTIONAL DIFFERENCES AMONG AQUIFICALES LINEAGES IN YNP Fixation of carbon dioxide**

A detailed inventory of genes coding for carbon dioxide fixation and various oxidation/reduction pathways indicated major functional differences among the three primary Aquificales lineages observed in this study (**Table 2**). The reverse-TCA pathway is thought to be the earliest CO2 fixation process used by microorganisms (Wachtershauser, 1990; Hugler et al., 2005). Although members of the Aquificales have been shown to utilize this pathway (Aoshima et al., 2004; Ferrera et al., 2007), the enzymes associated with the key step (citrate cleavage) differ among families of the cultured Aquificales. For example, *Sulfurihydrogenibium* spp. (Family

**FIGURE 6 | Hierarchical cluster analysis of relative gene abundances in theTIGRFAM role category "ElectronTransport" across six Aquificales streamer communities**. TIGRFAMs with low variation across the sites were removed before the clustering to retain ∼50 of the most variable families.

Subunits of the protein complexes were only represented by one representative TIGRFAM family. Pearson correlation was used as the distance measure for average linkage agglomerative clustering. Sites cluster consistent with pH and the presence or absence of sulfide.

 genes code for proteins with high specificity for possible

2No genes were found for nitrification (amoA), denitrification (e.g. nirK, nirS, norB), methanogenesis or methanotrophy (e.g. mcrA), sulfate reduction (e.g. dsrAB), or arsenate reduction (e.g. arrAB) within the Aquificales populations represented in these six sites.

3Includes Mo-pterin proteins similar to sreA and arrA.

 **a** Hydrogenothermaceae) catalyze citrate cleavage using ATP citrate lyase, which includes a large and small subunit (AclA and AclB, respectively). Conversely,*Thermocrinis* spp. (Family Aquificaceae) catalyze citrate cleavage using two separate enzymes, citryl-CoA synthetase (Ccs) and citryl-CoA lyase (Ccl). We detected the *aclB* gene in both *Sulfurihydrogenibium* sites (MHS\_10 and CS\_12), but not in any sites containing strictly Aquificaceae (DS\_9, OSP\_14, OS\_11 and BCH\_13; **Figure 7**). Although CS\_12 contained both *Thermocrinis* and *Sulfurihydrogenibium* populations, all *aclB* genes identified from this site grouped unambiguously with the *Sulfurihydrogenibium* sp. (**Figure 7**, bootstrap support = 92%).

Genes coding for citryl-CoA synthetase (*ccs*A) and citryl-CoA lyase (*ccl*) were found in both the *Hydrogenobaculum* (DS\_9 and OSP\_14) and *Thermocrinis*-like (OS\_11, CS\_12 and BCH\_13) populations (**Figures 8A,B**). Therefore, these members of the Aquificaceae appear to fix CO2 using this alternative citrate cleavage mechanism. Deduced protein sequences (CcsA and Ccl) from OS\_11, BCH\_13, DS\_9, and OSP\_14 grouped in distinct clades consistent with the pH differences among sites, as well as the different genera observed within this family (i.e., DS and OSP versus OS and BCH). The *ccs*A gene copy detected in CS\_12 grouped with similar entries in OS\_11 and BCH\_13, and was confirmed to come from the sub-dominant *Thermocrinis*-like population in this site. *Sulfurihydrogenibium*-like populations in MHS\_10 and CS\_12 lacked citryl-CoA synthetase (*ccs*A) genes, and instead contained one copy of a succinyl-CoA synthetase (**Figure 8A**). Differences in the gene neighborhood between the citryl-CoA and succinyl-CoA synthetase pathways of the Aquificaceae versus Hydrogenothermaceae suggest that the *Sulfurihydrogenibium* copy of succinyl-CoA synthetase from MHS is not involved in CO2 fixation. Annotation inconsistencies among these two fairly similar proteins (CcsA and succinyl-CoA synthetase) have made it difficult to make definitive assignments without visualizing the sequences in phylogenetic trees (**Figure 8A**) or other alignment tools, and these proteins are in fact thought to be related via gene duplication (Aoshima et al., 2004).

The r-TCA pathway has not actually been demonstrated in any cultured member of the *Hydrogenobaculum*, although field data on 14CO2 incorporation suggests that members of these communities fix CO2 at significant rates (Boyd et al., 2009). Genes found in DS\_9 and OSP\_14 are divergent relative to the *Thermocrinis* entries (**Figures 8A,B**), so it is unclear if these perform an identical function in both genera. However, other evidence that the r-TCA pathway is operative in the DS\_9 and OSP\_14 *Hydrogenobaculum*like populations includes two enzymes required for the reductive pathway: 2-oxoglutarate:ferredoxin oxidoreductase and fumarate reductase. Both genes are present in the *Hydrogenobaculum* sp. Y04AAS1 genome as well as the DS\_9 and OSP\_14 metagenomes.

*Oxidation of H2, reduced sulfur, and arsenite.* Other functional differences among the major Aquificales lineages include possible explanations for the specialization of these populations to specific geochemical environments. For example, *Hydrogenobaculum* populations from DS\_9 and OSP\_14 contain Group I Ni-Fe hydrogenases, but these genes are notably absent from the *Sulfurihydrogenibium* (MHS\_10 and CS\_12) and the *Thermocrinis* populations (OS\_11 and BCH\_13) (**Table 2**). The potential for H2

to serve as an electron donor for metabolism appears limited to the acidic sites where concentrations of aqueous H2 have been measured in the 50–100 nM range (Inskeep et al., 2005; Spear et al., 2005), while other sites contain lower H2 (**Table 1**). The *Hydrogenobaculum*-like populations were also the only Aquificales to contain genes coding for thiosulfate oxidase (*tqoAB*), often implicated in oxidation of thiosulfate in the order Sulfolobales (Friedrich et al., 2005). However, it is known that thiosulfate concentrations are considerably higher in circumneutral sulfidic sites, due to the greater stability of thiosulfate at intermediate pH (Xu et al., 1998; Nordstrom et al., 2005). The *Sulfurihydrogenibium-*like organisms may process thiosulfate through an abundance of rhodanese domain proteins known to be involved in sulfur-transferase reactions as well as SoxBC complexes (Friedrich et al., 2005). Moreover, all Aquificales lineages from each site contained either one or more copies of a highly conserved and syntenous gene complex thought to be important in the oxidation of reduced sulfur (i.e., sulfide and elemental S), and includes several hetero-disulfide reductases as well as other Fe-S proteins (*rhd*, *tusA*, *dsrE*, *hdrC*, *hdrB*, *hdrA*, *orf2*, *hdrC*, *hdrB)* (Ghosh and Dam, 2009). Each of the three lineages also contained *sqr* (sulfide:quinone reductase) genes (**Table 2**), which have been shown to code for proteins involved in the oxidation of dissolved sulfide to S0 or polysulfide chains, followed by electron transfer to the quinone pool through a flavin adenine dinucleotide (FAD) cofactor (Cherney et al., 2010). Even the *Thermocrinis* populations from OS\_11 and BCH\_13 exhibited potential for the oxidation of sulfide and elemental S, although it is unlikely that sufficient sulfide exists in these geothermal channels to support the growth of active "streamer" communities. The habitat range of *Thermocrinis*-like organisms in YNP includes high-pH (7–9) sulfidic channels (Inskeep et al., 2005; Hall et al., 2008; Planer-Friedrich et al., 2009) and thus may explain the presence of these genes in *Thermocrinis* assemblies. The HDR gene complexes appear highly conserved across numerous Aquificales and Sulfolobales, as well as acidophilic bacteria such as *Acidithiobacillus ferrooxidans* (Ghosh and Dam, 2009; Inskeep et al., 2013).

The oxidation of arsenite to arsenate is highly exergonic (ranging from 50 to 60 kJ/mol electron) in these geothermal systems (Inskeep et al., 2005) and has been shown to serve as a sole electron donor in several unrelated bacteria (D'Imperio et al., 2007; Santini et al., 2007). Consequently, it is interesting that both the *Hydrogenobaculum* and *Thermocrinis*-like organisms from sites DS\_9, OSP\_14, OS\_11, and BCH\_13 contain full copies of the arsenite oxidase Mo-pterin subunit I (*aio*A, also abbreviated *aro*A, *aso*A, *aox*B in prior work). *Thermocrinis* and *Hydrogenobaculum* spp. oxidize considerable amounts of arsenite in acidic to circumneutral springs in YNP (Macur et al., 2004; Inskeep et al., 2005; Hamamura et al., 2009), corresponding to the oxygenation of geothermal outflow channels and correlation with *aio*A expression in these same habitats (Clingenpeel et al., 2009; Hamamura et al., 2009). It is possible that members of the Aquificales gain energy from the oxidation of arsenite *in situ*. However, this has not been established in culture (Donahoe-Christiansen et al., 2004). All Aquificales populations from the current study also contained evidence for arsenic detoxification, including the potential to reduce arsenate via ArsC and extrude arsenite via an efflux pump (AsrB), as well as to methylate arsenite via methyl transferases (ArsM) (Bentley and Chasteen, 2002; Mukhopadhyay et al., 2002).

*Respiratory processes.* The presence of terminal oxidase complexes in each of the Aquificales populations suggests that these organisms all respire oxygen. However, the distribution of different types of subunit I heme Cu oxidases (HCOs) across the Aquificales lineages (**Figure 9**) is not consistent with the distribution of CO2 fixation genes. This observation invokes a different evolutionary history of carbon dioxide fixation versus aerobic respiration among these Aquificales lineages. The *Sulfurihydrogenibium* and *Hydrogenobaculum-*like organisms contain similar Type C-Cbb3 HCOs (Garcia-Horsman et al., 1994) despite the phylogenetic distance between these organisms (i.e., different families). Conversely, the *Thermocrinis* spp. from OS\_11 and BCH\_13 contained multiple copies (at least two distinct copies per site) of Type A HCOs (Pereira et al., 2001), which are phylogenetically related to HCOs from the majority of aerobic bacteria.

*Hydrogenobaculum* (DS\_9,OSP\_14) and *Sulfurihydrogenibium*like (MHS\_10, CS\_12) organisms were the dominant Aquificales in sites containing high levels of dissolved sulfide (low dissolved O2). Consequently, the Cbb3 cytochromes associated with these populations likely bind O2 with greater efficiency compared to Type A HCOs of the *Thermocrinis*-like populations, and is consistent with known properties of these proteins (Garcia-Horsman et al., 1994). The Cbb3-type heme-copper oxidases are found only in several groups of bacteria, especially the Aquificales and Proteobacteria (312 of 365 sequences; Sousa et al., 2011), and appear to represent the evolution of a separate respiratory complex in low O2 environments (Garcia-Horsman et al., 1994). In addition, the *Hydrogenobaculum* and *Sulfurihydrogenibium*-like organisms have genes coding for a *bd*-ubiquinol oxidase (*cydAB*), also thought to function under lower O2 concentrations (Jünemann, 1997). The *Thermocrinis*-like organisms in OS\_11 and BCH\_13 show no evidence of the *bd*-ubiquinol oxidases, consistent with the higher oxygen levels and low sulfide in these sites.

A detailed survey of other respiratory processes (dissimilatory reduction) suggests that the reduction of elemental sulfur and/or polysulfide is important in some members of the Aquificales. Although some sequenced Aquificales isolates contain a nitric oxide reductase (NorB) (one of the steps required for complete denitrification from nitrate to N2), *nar*G, *nor*B, *nir*K, or *nir*S like sequences were not found in the YNP Aquificales. Moreover, no evidence was found for dissimilatory reduction of sulfate or sulfite (*dsr*AB), arsenate (*arr*AB), or CO2 (*mcr*A, methanogenesis) in any of the Aquificales lineages. However, the organisms found in sulfidic channels (i.e., *Hydrogenobaculum* and *Sulfurihydrogenibium*) all contain sulfur reductases (*sre*A) or polysulfide reductase (*psr*A), as well as tetrathionate reductases (*ttr*A, another DMSO-Mopterin). Consequently, the Aquificales lineages that inhabit sulfidic channels under low O2 tensions may require electron shuttling to reduced sulfur instead of, or in addition to, O2. Genes for sulfur reduction were notably absent in the *Thermocrinis-*like populations of OS\_11 and BCH\_13, which correlates with the lack of elemental sulfur in these springs.

#### **DISCUSSION**

Phylogenetic and functional analysis of metagenome sequence from three major types of high-velocity, filamentous "streamer" communities revealed three lineages of Aquificales, whose metabolic potential correlated primarily with pH and sulfide and/or elemental sulfur. Sites with low-pH (pH 3–3.5) and high-sulfide contained *Hydrogenobaculum* spp., whereas higher-pH sites were dominated by either *Sulfurihydrogenibium* spp. (high-sulfide) or *Thermocrinis*-like (low sulfide) populations. *Calcite Springs* (CS\_12) also hosted a minority *Thermocrinis*-like population and was the only site here to contain two major Aquificales genera. This is consistent with previous 16S rRNA gene diversity surveys that have generally found only minor overlap in the distribution of different Aquificales across YNP geothermal environments (Reysenbach et al., 2005; Hall et al., 2008; Hamamura et al., 2009). *Thermocrinis* organisms have also been observed in sulfidic channels at higher-pH, near 9 (Planer-Friedrich et al., 2009).

Metagenome sequence assemblies for each of the three major Aquificales lineages resulted in total scaffold sizes that approach full genomes, and which represent "consensus sequence" or "pan-genomes" of these populations (Medini et al., 2005). Sequence variability of highly related populations within a field site may contribute to incomplete assembly, and "closure" of these *de novo* assemblies would require considerable manual effort, as well as additional sequencing to close gaps. Sequence heterogeneity within individual Aquificales populations was observed using AMPHORA to detect the number of single-copy genes. Although the streamer community from MHS\_10 was dominated by what appears to be a fairly homogeneous population type of *Sulfurihydrogenibium* sp., other sites exhibited greater variability within the primary Aquificales population. For example, the *Thermocrinis*like populations from *Octopus*, *Bechler*, and/or *Calcite Springs* all revealed higher numbers of numerous single-copy genes, which suggests greater heterogeneity of these populations *in situ*. Additional sequence coverage of these populations may result in less single-copy gene variability, and future efforts will be necessary to clarify sources of this variability. Importantly, the sequence assemblies generated in this study provide a foundation for future efforts to determine the types and rates of genetic change in these same sites.

In addition to the abundant Aquificales populations, *de novo* assemblies were also obtained for several novel bacterial and archaeal lineages, although at lower coverage. These cocommunity members provide an interesting comparative study in

**FIGURE 9 | Protein tree of heme-copper oxidases (subunit I of terminal oxidase complex)**. Metagenome entries are highlighted by site, and labels correspond to the major types of heme-copper oxidases observed in different Aquificales "streamer" communities. Similar heme-copper oxidases are found in sites OS\_11 and

BCH\_13, and these are significantly different than the terminal oxidases found in Hydrogenobaculum and or Sulfurihydrogenibium from sites DS\_9, OSP\_14, MHS\_10 and CS\_12 (neighbor-joining tree constructed using nitric oxide reductase (NorB) as the out group; 1000 bootstraps).

their own right given that each Aquificales community exhibited a different assemblage of interacting community members. For example, the acidic sites contained different archaeal populations representing specialization on Fe(II) (e.g., *Metallosphaera*-like) versus reduced sulfur (Thaumarchaeota- and Thermoplasmataleslike), while the higher-pH (pH ∼8) sites contained several novel bacterial lineages, and a population of Thermoproteaceae (G + C ∼ 60%) in OS\_11 and BCH\_13. Some of the novel bacterial populations present in OS\_11 and BCH\_13 were related to unclassified 16S rRNA gene sequences described previously (Reysenbach et al., 1994; Blank et al., 2002; Hall et al., 2008). Metagenomics provides a promising opportunity to gain insight into the metabolic potential of these novel populations, although additional sequencing will be required to bring the coverage up to levels suitable for building reasonable *de novo* assemblies (e.g., >2× coverage).

Global protein (TIGRFAM) analysis coupled with PCA and hierarchical clustering of all Aquificales "streamer" communities demonstrated important linkages among geochemistry, the presence of distinct phylotypes, and their metabolic genes. Clustering of TIGRFAMs specific to ET demonstrated specific linkages between individual phylotypes and respiratory pathways that are consistent with the strong influence of geochemistry on community structure (especially pH and sulfide). Perhaps one of the more interesting findings in the study is the degree to which co-community members varied across these high-velocity streamchannel habitats. Even site pairs containing the same Aquificales phylotype contained considerably different interacting populations. Factors controlling the variation in community structure in different sites containing the same Aquificales phylotype can also be shown to track with geochemical conditions. Such is the case when comparing the hypoxic elemental sulfur habitats (DS\_9) to the more oxic Fe(III)-oxide (OSP\_14) streamer communities, samples that both contained a very similar *Hydrogenobaculum* population. The archaeal co-community members in DS\_9 are likely anaerobic (or microaerophilic) populations compared to the aerobic Fe(II)-oxidizing Sulfolobales in OSP\_14.

Detailed analysis of functional genes present in the three Aquificales lineages revealed several examples of divergence that were likely driven by environmental selection and lateral transfer events. These events have resulted in inconsistent patterns in phylogeny between highly conserved genes and those coding for various functional processes. Characteristics of the r-TCA cycle make it a good candidate for analysis of early steps important in the evolution of autotrophy. The pathway's auto-catalytic nature and role in central C metabolism may reflect its importance in the evolution of associated anabolic and oxidative pathways (Hugler et al., 2005). The key ATP-dependent step of acetyl CoA synthesis is catalyzed by different proteins in the Hydrogenothermaceae versus Aquificaceae (AclB versus CcsA and Ccl, respectively; Aoshima et al., 2004), suggesting different evolutionary histories of these lineages with respect to CO2 fixation. The only Aquificales in this study to contain *acl*B were the *Sulfurihydrogenibium*-like organisms (Family Hydrogenothermaceae) present in MHS and CS (**Figure 7**). The *Hydrogenobaculum* and *Thermocrinis* populations (Family Aquificaceae) in the other sites contain *ccs*A and *ccl* genes rather than *acl*B (**Figures 8A,B**). The r-TCA pathway catalyzed by Ccs and Ccl is believed to be the ancestral pathway among the Aquificales

(Hugler et al., 2005), however, our understanding of this pathway is not complete as evinced by inconsistent annotation of both the citryl-CoA synthetase (*ccs*A) and the citryl-CoA lyase (*ccl*).

Although the Aquificales occupy diverse geochemical habitats, they generally flourish in zones of shallow, high-velocity, turbulent spring water where dynamic mixing occurs to create disequilibria between hypoxic and oxic conditions. The exchange of atmospheric oxygen and subsequent reaction with dissolved sulfide present in thermal waters is one of the more important geochemical processes occurring within sulfidic geothermal channels (Inskeep et al., 2005; Nordstrom et al., 2005). The disequilibrium between oxygen and sulfide establishes conditions suitable for microaerophiles growing on sulfide, elemental sulfur or thiosulfate (Reysenbach et al., 2009). Metagenome sequence of the YNP Aquificales clearly indicates the potential for the oxidation of reduced sulfur species using a variety of S-oxidation pathways coupled with Type C-Cbb3 or *bd*-ubiquinol terminal oxidase complexes, especially in the *Hydrogenobaculum* and *Sulfurihydrogenibium-*like organisms detected in sulfidic systems. The *Thermocrinis* organisms present in OS\_11 and BCH\_13 contain Type A-HCOs (**Figure 9**), indicating their functional divergence from the other Aquificales genera with respect to oxygen. There is evidence that some Aquificales have copies of both types of HCOs as is noted for *Hydrogenivirga* spp. Consequently, it is possible that subsequent evolution in specific habitat types separated lineages containing only the Type C or the Type A HCO. Interestingly, the *Hydrogenobaculum* and *Sulfurihydrogenibium* populations both contain the cbb3-Type C HCOs even though they are members of different families. The *Hydrogenobaculum* and *Thermocrinis* organisms are from the same family, but do not share the same respiratory complexes, although they do share similarity in r-TCA proteins important in CO2 fixation. The major difference in metabolic processing of CO2 and O2 between these lineages provides an excellent opportunity for relating their evolutionary histories to paleobiological events and the timing of radiation relative to the "Great Oxidation Event" (Canfield, 2005; Anbar et al., 2007; Konhauser, 2009).