Expanded Diversity and Phylogeny of mer Genes Broadens Mercury Resistance Paradigms and Reveals an Origin for MerA Among Thermophilic Archaea.

Mercury (Hg) is a highly toxic element due to its high affinity for protein sulfhydryl groups, which upon binding, can destabilize protein structure and decrease enzyme activity. Prokaryotes have evolved enzymatic mechanisms to detoxify inorganic Hg and organic Hg (e.g., MeHg) through the activities of mercuric reductase (MerA) and organomercury lyase (MerB), respectively. Here, the taxonomic distribution and evolution of MerAB was examined in 84,032 archaeal and bacterial genomes, metagenome assembled genomes, and single-cell genomes. Homologs of MerA and MerB were identified in 7.8 and 2.1% percent of genomes, respectively. MerA was identified in the genomes of 10 archaeal and 28 bacterial phyla previously unknown to code for this functionality. Likewise, MerB was identified in 2 archaeal and 11 bacterial phyla previously unknown to encode this functionality. Surprisingly, homologs of MerB were identified in a number of genomes (∼50% of all MerB-encoding genomes) that did not encode MerA, suggesting alternative mechanisms to detoxify Hg(II) once it is generated in the cytoplasm. Phylogenetic reconstruction of MerA place its origin in thermophilic Thermoprotei (Crenarchaeota), consistent with high levels of Hg(II) in geothermal environments, the natural habitat of this archaeal class. MerB appears to have been recruited to the mer operon relatively recently and likely among a mesophilic ancestor of Euryarchaeota and Thaumarchaeota. This is consistent with the functional dependence of MerB on MerA and the widespread distribution of mesophilic microorganisms that methylate Hg(II) at lower temperature. Collectively, these results expand the taxonomic and ecological distribution of mer-encoded functionalities, and suggest that selection for Hg(II) and MeHg detoxification is dependent not only on the availability and type of mercury compounds in the environment but also the physiological potential of the microbes who inhabit these environments. The expanded diversity and environmental distribution of MerAB identify new targets to prioritize for future research.

Introduction

Mercury (Hg) is one of the most toxic heavy metals and is a global pollutant emitted by volcanoes, forest fires, and anthropogenic activities such as fossil fuel combustion and smelting and oxidative leaching of pyritic ores (Rooney, 2007; Obrist et al., 2018). The global distribution of Hg poses severe consequences to human and environmental health (Clarkson and Magos, 2006), and the Hg load to environments is only increasing due to increased fossil fuel emissions, changes in land use, and artisanal gold mining (Obrist et al., 2018; Zolnikov and Ramirez Ortiz, 2018). Volatile mercury [Hg(0)] is the primary form of Hg emitted to the atmosphere where it can then be oxidized to Hg(II) and returned to terrestrial and aquatic ecosystems through wet deposition (Selin et al., 2007). Subsequently, Hg(II) can be buried in sediments, reduced to Hg(0) and re-emitted to the atmosphere, or methylated to form organic methylmercury (MeHg) (Selin et al., 2007; Grégoire and Poulain, 2018). MeHg is a highly potent neurotoxin that can enter food webs (Gilmour and Henry, 1991; Watras et al., 1995; Benoit et al., 1998) through uptake by unicellular microbes (Mason et al., 1996; Pickhardt and Fisher, 2007). MeHg can bioaccumulate and its concentration biomagnify in aquatic (Watras et al., 1998; Boyd et al., 2009) and terrestrial (Rimmer et al., 2010) food chains, ultimately impacting the health of animals and plants (Sunderland, 2007; Zhang et al., 2010).

Microorganisms are the primary driver of Hg methylation through the activity of a corrinoid protein, HgcA, and a ferredoxin, HgcB (Gilmour et al., 2013; Parks et al., 2013). Microbial Hg methylation occurs mainly under anaerobic conditions through the activities of iron-reducing bacteria (Fleming et al., 2006), sulfate-reducing bacteria (Compeau and Bartha, 1985; Gilmour et al., 1992), methanogenic archaea (Hamelin et al., 2011; Gilmour et al., 2018), and fermentative bacteria (Parks et al., 2013). Recent bioinformatic investigations of the distribution of hgcAB genes in metagenomes (Podar et al., 2015) and in isolate genomes and metagenome assembled genomes (McDaniel et al., 2020) showed that microbial Hg methylation potential is encoded in the genomes of a variety of organisms inhabiting anoxic environments with temperatures < 70°C but was less abundant in environments with temperatures > 70°C and notably absent from human and mammalian microbiomes.

Functionalities encoded by the mer operon are the primary mechanisms microbial cells use for Hg resistance and detoxification (Silver and Hobman, 2007; Lin et al., 2011). At its core, the mer operon encodes a homodimeric flavin-dependent disulfide oxidoreductase, termed mercuric reductase (MerA), that functions to reduce Hg(II) to volatile Hg(0) that can then diffuse out of the cell (Boyd and Barkay, 2012; Figure 1). The operon may also code for organomercury lyase (MerB), that catalyzes the protonolytic cleavage of the C-Hg bond in organomercury compounds, among them MeHg (Boyd and Barkay, 2012). This reaction yields a reduced organic moiety, that in the case of MeHg is methane, and Hg(II); the latter is then reduced to Hg(0) through the activity of MerA (Boyd and Barkay, 2012). In addition to MerAB, mer operons may code for the periplasmic protein MerP, inner membrane spanning proteins MerT, MerC, MerE, MerF, and MerG that transport Hg(II) to the cytoplasmic MerA (reviewed in Boyd and Barkay, 2012). mer operons may also code for transcriptional regulators such as the repressor/activator MerR or anti-activator MerD (reviewed in Barkay et al., 2003; Lin et al., 2011). Through the combined demethylation of MeHg by MerB and decreased concentration of Hg(II) in the cytoplasm by its reduction to gaseous Hg(0) via MerA, the mer operon effectively partitions Hg to the gaseous phase, allowing for microbial growth (Barkay et al., 2003; Boyd and Barkay, 2012).

FIGURE 1

The mer detoxification system. (A) A generic mer operon. (B) The cellular mer-encoded mercury detoxification mechanisms. The outer cell wall is depicted by a broken line illustrating that not all microbes have an outer membrane; broken line arrows depict diffusion; solid line arrows indicate transport or transformations; L = ligand with subscripts denoting the ligand type. The colors of various Mer proteins correspond with the colors of the genes that encode these proteins in panel A. Figure adapted from Lin et al. (2011); Boyd and Barkay (2012).

In 2012, a database of publicly available genomes was compiled that included 272 mer operons from the genomes of a diverse array of 246 bacterial and archaeal isolates (Boyd and Barkay, 2012). Phylogenetic analyses of MerA homologs compiled in that study, as well as in a previous study (Barkay et al., 2010), found that the earliest branching MerA were from thermophilic bacteria (Aquificales) and that merA was acquired in Archaea through lateral gene transfer. Further, it was shown that the complexity of mer operons, as indicated by the presence of other mer genes found linked to MerA in an operonic structure, increased with evolutionary time through recruitment of additional functionalities. As a result, the efficacy and function of the mer detoxification machine was improved including a broadened spectrum of mercurials that could be detoxified by the addition of MerB (Boyd and Barkay, 2012). Targeted studies of the diversity and distribution of MerA-encoding genes have since shown a greater diversity in natural environments than was found in the complete genomes of isolates. This includes diverse merA detected in studies of high Arctic snow, freshwater, and sea-ice brine (Møller et al., 2014), in Antarctic sea ice (Gionfriddo et al., 2016), in deep and hot brine environments of the Red Sea (Sayed et al., 2014; Maged et al., 2019), in bioreactors (Shi et al., 2019), and in mercury-contaminated rice paddies in China (Vishnivetskaya et al., 2018), among others. While merA is estimated to be roughly ten-times more abundant than merB in genomes of isolates (Boyd and Barkay, 2012), more recent metagenomics studies have shown that the number of genomes that encode merB is roughly equivalent to those that encode hgcAB (Christensen et al., 2019). This is a significant finding considering the global prevalence and distribution of organisms with the ability to generate MeHg (despite it being a relatively rare microbial trait) through the activity of HgcAB (Podar et al., 2015; McDaniel et al., 2020) and the potential for MerB-encoding organisms to balance MeHg production in the environment through demethylation.

Since the original mer database was compiled in 2012 (Boyd and Barkay, 2012), significant advances in sequencing technologies, including the broad application of metagenomics and single cell genomics in diverse environmental systems (Bowers et al., 2017; Parks et al., 2017, 2020), have substantively increased the known diversity of microbial life (Brown et al., 2015; Chen L.X. et al., 2018). Yet, the diversity of MerAB protein homologs in these newly discovered microbial genomes remains uncharacterized. In the present study, bioinformatic approaches were applied to publicly available isolate genomes, metagenomic assembled genomes (MAGs), and single cell genomes (SCGs). Genomes, MAGs and SCGs were obtained from the Department of Energy (DOE) Integrated Microbial Genomes/Microbiomes (IMG/M) database (as of May 2020) to investigate the taxonomic distribution and diversity of MerAB, the centerpieces of the mer Hg detoxification and MeHg degradation machinery (Boyd and Barkay, 2012). Phylogenetic approaches were applied to reconstruct the evolutionary history of MerA homologs and to provide insight into when MerB was first recruited to the mer operon. Finally, potential alternative mechanisms allowing for Hg(II) detoxification in cells that were found to encode MerB but not MerA are proposed to explain how toxic effects of Hg(II) generation from demethylation of MeHg might be overcome.

Materials and Methods

Generation of a MerA and MerB Homolog Database

Mercuric reductase (MerA) protein homologs were compiled from publicly available genomes deposited in the DOE IMG/M database (Chen I.M.A. et al., 2018) including isolate genomes, MAGs, and SCGs, in May 2020, using BLASTp within IMG (Markowitz et al., 2006) and MerA from the bacterial transposon Tn501 [CAA77323.1; Stanisich et al., 1977] and the archaeal taxon Saccharolobus solfataricus P2 [AAK42805.1; Schelert et al., 2004] as queries. Additional sequences were identified with the Gene Search function within IMG/M with input query terms “mercury(II) reductase,” “MerA,” and “mercuric reductase.” All putative MerA sequences were aligned with the PROMALS3D structural alignment tool (Pei et al., 2008) using the protein crystal structure of Tn501 MerA [PDB ID: 1ZK7 (Ledwidge et al., 2005)] as the alignment thread.

Aligned putative MerA were manually examined for the presence of sequence signatures that have been experimentally shown to be essential for MerA activity, as previously described (Barkay et al., 2003; Boyd and Barkay, 2012). Specifically, homologs were screened for the presence of the vicinal cysteine pair at the carboxy terminus (aligned positions 628 and 629; numbering in reference to MerA from Bacillus cereus RC607 [BAB62433.1; Wang et al., 1989; Gupta et al., 1999] and for the conserved cysteine pair at positions 207 and 212 in the redox active site, tyrosine at position 264 (Rennex et al., 1993), and tyrosine at position 605 for bacterial MerA and phenylalanine at position 605 for archaeal MerA (Simbahan et al., 2005).

The IMG/M database was also queried for the presence of alkylmercury or organomercurial lyase (MerB) using BLASTp and MerB from the Escherichia coli plasmid R831b [AAB49639.1; Pitts and Summers, 2002], Halobacterium sp. DL1 [AHG04994.1] and Candidatus Nitrososphaera evergladensis SR1 [AIF84719.1; Zhalnina et al., 2014] as queries. Additional MerB homologs were identified with the Gene Search function within IMG/M with input query terms “merB,” “alkylmercury lyase,” “alkylmercury lyase-like,” “organomercury lyase,” and “organomercurial lyase.” Homologs were aligned with the PROMALS3D structural alignment tool specifying the protein crystal structure of MerB from plasmid R831b [PDB ID: 3F0P (Lafrance-Vanasse et al., 2009) as the alignment thread. Aligned MerB homologs were examined for the presence of the following sequence signatures: cysteine at position 96 (numbering in reference to MerB from R831b), aspartic acid at position 99, and cysteines at position 159 and 117 (Pitts and Summers, 2002; Lafrance-Vanasse et al., 2009). Pairwise sequence similarity among MerA or MerB homologs encoded in the same genome were determined using the Sequence Manipulation Suite: Version 2 (Stothard, 2000).

Phylogenetic Analysis of MerA and MerB

MerA protein homologs were aligned using Clustal Omega version 1.2.4 (Sievers et al., 2011), using default settings with dihydrolipoamide dehydrogenase from Magnetospirillum magneticum AMB-1 (YP_423326), Thermus thermophilus HB27 (YP_005669), Pseudomonas fluorescens Pf0-1 (WP_011336663.1), Desulfurococcus amylolyticus DSM 16532 (AFL66567.1), and Saccharolobus solfataricus (SAI85240.1) serving as outgroups. N terminal “NmerA” sequences were trimmed from the alignment block as previously described (Barkay et al., 2010; Boyd and Barkay, 2012). The aligned homologs were subjected to Modelfinder included in IQ-TREE vers. 1.6.12 (Kalyaanamoorthy et al., 2017) to identify the optimal substitution model and parameters for the phylogenetic reconstruction of MerA. The phylogeny of MerA was evaluated with IQ-TREE multicore version 1.6.11 (Nguyen et al., 2015; Hoang et al., 2018) with the LG amino acid substitution matrix, empirically counted frequencies from the alignment, a discrete four category gamma substitution model (gamma shape parameter = 0.9583), and a defined proportion of invariant sites of 0.0080, specifying 1000 ultrafast bootstraps. The phylogenetic tree was visualized in ITOL (Letunic and Bork, 2019). The MerA phylogenetic tree is available in newick format in Supplementary Dataset 1. All aligned MerA homolog sequences used for the construction of the phylogenetic tree, along with the plasmid MerA homologs not used for the phylogenetic tree construction, are available in fasta format in Supplementary Dataset 2.

Compiled MerB protein homologs were aligned with Clustal Omega version 1.2.4 (Sievers et al., 2011), specifying default settings. The aligned homologs were subjected to Modelfinder included in IQ-TREE version 1.6.12 as mentioned above. The phylogeny of MerB was evaluated with IQ-TREE multicore version 1.6.11 as mentioned above but with gamma shape parameter = 1.6454, defined proportion of invariant sites of 0.0000, and 1000 ultrafast bootstraps. The phylogenetic tree was visualized in ITOL. The MerB phylogenetic tree is available in newick format in Supplementary Dataset 3. All aligned MerB homolog sequences are available in fasta format in Supplementary Dataset 4.

Results and Discussion

Overview of the Genomic Database Queried for MerA and MerB

A total of 84,032 isolate genomes, MAGs, and SCGs (heretofore referred to as genomes) available in the IMG/M database as of May 2020 were queried for homologs of MerA and MerB. Of the 84,032 genomes queried, 80,925 were bacterial, 1959 were archaeal, and 1148 were plasmid genomes. Furthermore, 70,633 were isolate genomes (including plasmid genomes), 9590 were MAGs, and 3809 were SCGs. MerA and/or MerB homologs were identified in 7161 isolate genomes, 480 MAGs and 76 SCGs (Supplementary Table 1).

Taxonomic Distribution of MerA Encoding Genes

A total of 7544 MerA homologs were identified (Supplementary Table 1) and these were distributed among 6574 archaeal and bacterial genomes (Table 1; Supplementary Table 1), with 813 of these genomes encoding for multiple (up to 8) homologs (Table 2). The 7544 MerA homologs were derived from 6063 isolate and plasmid genomes, 446 MAGs, and 65 SCGs (Supplementary Table 1).

TABLE 1

Number of genomes that encode MerA or MerB protein homologs and their distribution at the Phylum level in the database constructed in 2012 (Boyd and Barkay, 2012) and here in 2020.

	Number of genomes encoding MerA homologs		Number of genomes encoding MerB homologs*
Phylum/Class	2020	2012	2020	2012
Archaea	271	24	11	−
Candidatus Bathyarchaeota	1	−	−	−
Candidatus Diapherotrites	1	−	−	−
Candidatus Geoarchaeota	4	−	−	−
Candidatus Geothermarchaeota	2	−	−	−
Candidatus Heimdallarchaeota	1	−	3	−
Candidatus Micrarchaeota	22	−	−	−
Candidatus Odinarchaeota	1	−	−	−
Candidatus Parvarchaeota	10	−	−	−
Crenarchaeota	139	15	−	−
Euryarchaeota	65	9	8	−
Nanoarchaeota	1	−	−	−
Thaumarchaeota	10	−	−	−
unclassified	14	−	−	−
Bacteria	6274	253	1778	40
Acidobacteria	7	−	−	−
Actinobacteria	958	36	468	13
Aquificae	7	3	−	−
Armatimonadetes	1	−	−	−
Bacteroidetes	131	5	3	4
Balneolaeota	4	−	−	−
Candidate division GAL15	1	−	−	−
Candidate division NC10	−	−	1	−
Candidate division Zixibacteria	−	−	1	−
Candidatus Atribacteria	1	−	−	−
Candidatus Azambacteria	1	−	−	−
Candidatus Bipolaricaulota	1	−	−	−
Candidatus Blackburnbacteria	1	−	−	−
Candidatus Daviesbacteria	7	−	−	−
Candidatus Giovannonibacteria	10	−	−	−
Candidatus Gottesmanbacteria	1	−	−	−
Candidatus Kaiserbacteria	1	−	−	−
Candidatus Kryptonia	1	−	−	−
Candidatus Microgenomates	2	−	−	−
Candidatus Pacebacteria	1	−	−	−
Candidatus Parcubacteria	2	−	−	−
Candidatus Peregrinibacteria	5	−	−	−
Candidatus Roizmanbacteria	2	−	−	−
Candidatus Rokubacteria	2	−	1	−
Candidatus Sumerlaeota	1	−	−	−
Chlamydiae	1	−	1	−
Chloroflexi	32	1	7	−
Cyanobacteria	1	−	−	−
Deferribacteres	2	−	−	−
Deinococcus-Thermus	14	2	2	−
Firmicutes	960	42	797	5
Fusobacteria	−	−	1	−
Gemmatimonadetes	−	−	1	−
Ignavibacteriae	15	−	−	−
Nitrospirae	18	1	6	−
Planctomycetes	1	−	−	−
Proteobacteria	4050	161	478	18
Alphaproteobacteria	668	25	84	3
Betaproteobacteria	412	24	72	2
Gammaproteobacteria	2869	104	298	13
Deltaproteobacteria	37	8	10	−
Epsilonproteobacteria	1	−	−	−
Zetaproteobacteria	27	−	4	−
Acidithiobacillia	25	−	6	−
Hydrogenophilalia	2	−	1	−
Oligoflexia	6	−	2	−
Candidatus Muproteobacteria	1	−	−	−
Unclassified	2	−	1	−
Rhodothermaeota	1	−	−	−
Spirochaetes	3	−	5	−
Synergistetes	−	−	3	−
Tenericutes	11	1	−	−
Verrucomicrobia	7	1	−	−
unclassified	10	−	3	−
Plasmid:Bacteria	29	−	6	−
Alphaproteobacteria	1	−	1	−
Betaproteobacteria	1	−		−
Gammaproteobacteria	25	−	4	−
Unclassified Proteobacteria	2	−	1	−
Ynclassified	−	3	−	3

TABLE 2

Abundance of archaeal, bacterial, and plasmid genomes that encode one or more MerA or MerB protein homologs.

MerA copies per genome	1	2	3	4	5	6	8	Total
Archaea	270	−	1	−	−	−	−	271
Bacteria	5497	668	74	29	2	3	1	6274
Plasmid:Bacteria	28	1	−	−	−	−	−	29
MerB copies per genome	1	2	3	4	8	Total
Archaea	11	−	−	−	−	11
Bacteria	1685	64	15	13	1	1778
Plasmid:Bacteria	5	1	−	−	−	6

Among the 1959 archaeal genomes examined, a total of 271 (13.8%) encoded MerA homologs and these were distributed among 12 phyla (Table 1). The 307 archaeal MerA protein homologs in the current database represents a roughly ten-fold increase over the number identified in the 2012 database (Table 1; Boyd and Barkay, 2012), and include homologs from an additional 10 phyla and candidate phyla. These include genomes affiliated with the Thaumarchaeota, Nanoarchaeota, and various candidate phyla including Candidatus Micrarchaeota, Candidatus Parvarchaeota, and Candidatus Geoarchaeota. While this is an expansion of our previous inventories, the added archaeal genomes are from taxa that generally inhabit hot springs, in particular acidic hot springs that are acidified through microbially-mediated O2-dependent oxidation of sulfur compounds (Colman et al., 2018). The presence of MerA encoding genomes from taxa in these environments is thus consistent with the aerobic conditions that have been previously suggested to enhance Hg toxicity (Barkay et al., 2010). Nearly all (99.6% of total) MerA encoding archaeal genomes coded for a single copy (Table 2) with the exception of Cuniculiplasma divulgatum S5 (Euryarchaeota; Diaforarchaea group) that encodes 3 copies. Two of the three copies of MerA encoded in C. divulgatum are highly similar (89.2% sequence identities), suggesting a possible gene amplification within this lineage. C. divulgatum S5 is an acidophilic, mesophilic, facultatively anaerobic heterotroph isolated from a sulfide ore mine (Golyshina et al., 2016). Sulfide ores often have trace Hg(II) (Monecke et al., 2016; Fallon et al., 2017) that, when liberated through oxidative weathering, might be expected to select for duplicons encoding for MerA to increase Hg(II) detoxification capabilities through a gene dosage effect (Kondrashov, 2012).

Among the 84,032 bacterial genomes examined, a total of 6274 (7.75%) encoded MerA homologs (Table 1). An additional 29 proteobacterial plasmid genomes encoded MerA homologs (Table 1), out of the 1148 plasmid genomes queried. The 7207 MerA homologs were distributed among 37 phyla, candidate phyla, or candidate divisions. The number of genomes that encode MerA in the present genomic database is roughly 25-fold greater than those identified in 2012 (Table 1; Boyd and Barkay, 2012) and includes MerA from 26 new bacterial phyla including Acidobacteria, Balneolaeota, and Candidatus Giovannonibacteria, among others (Table 1). Further, homologs of MerA were identified among several new proteobacterial classes such as the Epsilonproteobacteria (e.g., Campylobacter concisus), the iron-oxidizing Zetaproteobacteria, and the thermophilic Hydrogenophilalia. Of the 6274 bacterial genomes that encode MerA, 777 encode multiple copies (Table 2), including one genome (Cupriavidus metallidurans H1130) that encodes 8 copies. While C. metallidurans strain H1130 was isolated from the blood of a patient suffering from nosocomial septicemia (Monsieurs et al., 2013), other C. metallidurans strains are frequently isolated from environments impacted by mining, metallurgic, and chemical industries (Diels and Mergeay, 1990; Goris et al., 2001). All 8 copies of MerA in this genome are identical (100% sequence identities) suggesting recent duplication events within this lineage, perhaps a response to increase detoxification efficiency due to Hg(II) exposure.

Several MerA homologs contained residues that vary from characterized sequences, in particular, at position 605 (numbering in reference to MerA of Bacillus cereus RC607). Variants at this position were also detected in a previous characterization of MerA diversity (Boyd and Barkay, 2012), with tyrosine identified in bacterial MerA homologs and phenylalanine in archaeal homologs (Supplementary Figure 1). We detected 46 bacterial MerA homologs that feature a phenylalanine instead of a tyrosine at position 605 and 2 archaeal MerA homologs that feature a tyrosine instead of a phenylalanine in this same position (Supplementary Table 2). The 46 bacterial MerA variants are distributed among 15 phyla, while the 2 archaeal variants both belong to members of the Euryarchaeota (Supplementary Table 2). All but two of these MerA variants are the sole MerA homologs in their respective genomes, while 2 of these variants co-exist in the genome with a MerA homolog with the typical residue at that position (canonical MerA) (Supplementary Table 1). To date, no MerA with a noncanonical substitution in position 605 has been tested for activity. Demonstrating activity of such variants may lead to a better understanding of the role of this residue in Hg(II) reduction by MerA.

Taxonomic Distribution of MerB Encoding Genes

A total of 1936 alkylmercury lyase (MerB) homologs were identified in publicly available genomes of the IMG database. These 1936 MerB protein homologs were distributed among 1722 genomes, 61 MAGs, and 12 SCGs (Supplementary Table 1), with 94 genomes encoding multiple homologs (up to 8; Table 2). While several genomes analyzed herein were shown to encode multiple copies of MerB, to date no genome has been shown to encode multiple copies of HgcAB (Podar et al., 2015; McDaniel et al., 2020).

Among the 1959 archaeal genomes examined, 11 (0.6%) encoded MerB homologs and these were distributed across 2 phyla: Candidatus Heimdallarchaeota and Euryarchaeota (Table 1). Demethylation of MeHg has not yet been described among the Archaea and a previous study of available isolate genomes in 2012 did not identify MerB homologs among the MerA-linked loci in archaeal genomes (Boyd and Barkay, 2012). The detection of MerB homologs in archaeal genomes is consistent with the significant expansion of the taxonomic diversity of the mer operons in the present study. All archaeal MerB-encoding genomes encoded for a single MerB copy (Table 2). Interestingly, MerB homologs were not identified among genomes of members of the Crenarchaeota, the phylum that featured the largest number of MerA homologs (Table 1).

Among the 80,925 bacterial genomes examined, 1778 (2.4%) encoded MerB protein homologs and these were distributed across 15 different phyla, candidate phyla, and divisions (Table 1), including 11 new phyla and candidate divisions not reported to encode MerB in the 2012 genomic database (Boyd and Barkay, 2012). These include members of the Nitrospirae, Spirochetes, and Chloroflexi along with the proteobacterial classes Deltaproteobacteria, Zetaproteobacteria, and Acidithiobacillia, among others (Table 1). Additionally, 7 homologs were identified on 6 proteobacterial plasmid genomes, with 1 plasmid containing 2 homologs (Table 2). Studies on organomercury resistance and degradation of MeHg have largely been focused on proteobacterial (Pitts and Summers, 2002) and bacilli (Matsui et al., 2016) strains. The broad taxon distribution of MerB homologs reported here calls for an examination of enzyme activity in more organisms where this activity has not been tested before.

Among the 1795 MerB-encoding genomes, 165 encoded multiple copies, up to 8 in the genome of Bacillus cereus strain #17 (Supplementary Table 3; Table 2). This strain was isolated from a mouse gut microbiome sample (Böhm et al., 2015). These 8 homologs probably originated from 2 different ancestral genes, with each of these genes undergoing multiple (quadruple) duplication events and thus forming 2 groups that each comprise 4 identical (100% sequence identities) homologs; the 100% similarity within each group suggests that the duplication events took place recently. This strain further encodes one additional group of 4 identical copies of the MerB-like variant, 99Ser, which is described in detail below. Bacillus strains have been previously found to encode multiple merB genes, designated as merB1, merB2, merB3 (Gupta et al., 1999; Huang et al., 1999; Chien et al., 2010; Boyd and Barkay, 2012). Heterologous expression and biochemical characterization reveal Bacillus MerB1 to preferentially degrade alkylmercurials (methylmercury, ethylmercury, thimerosal, and p-chloromercuribenzoate) whereas Bacillus MerB3 preferentially degrades arylmercurials (p-chloromercuribenzoate) (Chien et al., 2010). It is likely that the multiple copies of MerB in the B. cereus strain #17 broadens the range of organomercurials that this strain can detoxifiy. The B. cereus #17 genome also contains 4 copies of merA. While little information is available on the mouse host from which this strain was isolated (Böhm et al., 2015), it is possible that the diet of this mouse was rich in organomercurials, leading to selection for a strain that encoded such a diverse array of MerAB. Alternatively, it is possible that this strain or a recent ancestor of this strain was exposed to environmental organomercurials that selected for this high level of redundant functionalities prior to its recruitment as a member of the mouse gut microbiome.

MerB protein homologs were initially identified and screened using specific amino acid signatures, as previously described (Pitts and Summers, 2002; Lafrance-Vanasse et al., 2009). This includes cysteine at position 96 (numbering in reference to MerB from plasmid R831b), aspartic acid at position 99, and cysteines at position 159 and 117. Apart from the MerB homologs identified, an additional number of MerB-like homologs were identified which differed in these sequence signatures. When a homolog differed in at least one of the aforementioned amino acid signatures, these homologs were designated as MerB-like variants 99Ser and 117alt (Supplementary Figure 1). In MerB-like variant 99Ser, serine is found at position 99 instead of aspartic acid, and in MerB-like variant 117alt, an alternative amino acid is present at position 117 instead of cysteine (Supplementary Figure 1; Supplementary Table 3). The 13 MerB-like 99Ser variants were identified in genomes that all contain at least one MerB with the aforementioned canonical sequence signatures; these genomes belonged to members of the bacterial phylum Firmicutes. Previous structural and biochemical characterization of representatives of the 99Ser MerB-like variants (Wahba et al., 2016) has shown that they have expanded range of affinities for metals relative to canonical MerB. Specifically, this variant has an affinity to bind Cu(II) that is displaced by Hg(II) when present. Hg(II) is a product of MerB-dependent demethylation of MeHg, and MerB binds Hg(II) so that it cannot diffuse into the cytoplasm; MerA then extracts the metal for reduction (Wahba et al., 2016). The physiological role of these variant enzyme with an affinity for other divalent metals is thus not well understood especially when a native MerB 99Ser, in B. megaterium MerB2, does not bind Cu(II) (Wahba et al., 2016). However, it is possible that they function as more generalized organic metalloid detoxification enzymes than canonical MerB. If so, strains carrying MerB-like 99Ser in addition to the canonical MerB may possess an expanded range of metal resistance.

The importance of the 117Cys to MerB function is not clear. Pitts and Summers (2002), who initially identified it as a conserved amino acid, attributed a structural rather than a catalytic role to this moiety (Barkay et al., 2003). Current models of MerB show 117Cys to be deeply buried in the hydrophobic core of the protein and removed from solvent, and thus interactions with putative substrates (Di Lello et al., 2004). Of the total 563 MerB-like 117alt variant homologs, 20 were identified in genomes with a canonical MerB and these genomes belong to Actinobacteria and Firmicutes (Supplementary Table 3). The remainder of the 563 MerB 117alt variants were the only MerB homologs detected in their corresponding genomes and they are distributed in Actinobacteria, Bacteroidetes, Chloroflexi, Deinococcus-Thermus, Firmicutes, Proteobacteria, and Spirochetes (Supplementary Table 3). Crystallographic, structural, and biochemical studies of MerB have focused on the essential role of Cys96, Asp99, and Cys159 residues (Lafrance-Vanasse et al., 2009; Wahba et al., 2016; Roy et al., 2020), likely due to their strict conservation and catalytic role. To the best of our knowledge, the activity and/or substrate range of variants with substitution at position 117 has not been tested; a C117S mutant was insoluble and could not be tested (Pitts and Summers, 2002). The widespread taxonomic distribution of this variant in genomes lacking a canonical MerB suggests that the activity of 117alt should be tested for a better understanding of the physiological and ecological role of this enzyme.

Co-distribution of MerA and MerB Homologs in Genomes

A total of 5685 genomes encoded mer operons that included MerA but not MerB. This suggests a decreased spectrum of Hg compounds that these cells can detoxify. Interestingly, a total of 906 genomes, slightly higher than 50% of all merB-carrying genomes (a total of 1795 genomes), encoded MerB protein homologs but no homolog of MerA (Figure 2; Supplementary Table 1). These included 10 draft genomes and 896 finished genomes indicating that this observation is unlikely to be an artifact of incomplete genome sequencing. Moreover, of the 906 MerB-encoding genomes that lacked a MerA homolog, 900 belonged to Bacteria, 1 was carried by a bacterial plasmid, and 5 belonged to Archaea. Of these 906 genomes, 870 are isolate genomes, 29 are MAGs, and 7 are SCGs. It should be noted that, in the case of MAGs, it is possible that a MerA-encoding plasmid could have been missed during the assembly process (Maguire et al., 2020). The presence of a plasmid-encoded MerA could be identified in future studies through application of specific algorithms that are designed to identify plasmid genomes in assembled metagenomic data (Rozov et al., 2017; Antipov et al., 2019). Of the 6574 MerA-encoding genomes, 249 encoded only 117alt MerB-like variants, 36 encoded 117alt MerB-like variants along with at least one canonical MerB, and 9 genomes encoded a MerB-like 99Ser variant along with a canonical MerB. Of the 895 genomes that encode a MerB homolog but lacked a MerA homolog, 732 (81.3%) belong to the genera Clostridioides, Mycobacteroides, Clostridium, Staphylococcus, Corynebacterium, Streptomyces, Paeniclostridium, Amycolatopsis, Thioalkalivibrio, Mycobacterium, Nocardia, or Enterobacter (Supplementary Table 1). These genera include aerobic, facultatively anaerobic, microaerophilic, and anaerobic representatives, indicating that the ability to demethylate MeHg via MerB without a dedicated strategy to reduce and volatilize Hg(II) via MerA is not necessarily related to the strains ability to use oxygen. Additionally, none of the genomes that encode only MerB code for homologs of HgcAB. Further, mapping of the presence of MerA on an unrooted MerB and MerB-like phylogeny failed to reveal a pattern in the distribution of organismal genomes that encode only MerB homologs versus those that encode both MerAB (Supplementary Figure 1). To the best of the authors’ knowledge, activities of organisms carrying only MerB or organisms that encode only MerB have yet to be reported (Boyd and Barkay, 2012; Grégoire and Poulain, 2018). Organomercury degradation without removal of Hg(II) may present challenges to the activity of the enzyme and the growth and survival of the organism. First, current models suggest a direct transfer of the Cys96-bound Hg(II) from MerB to MerA (Benison et al., 2004) raising questions about how MerB catalysis can continue in absence of MerA. Second, production of intracellular Hg(II) from MeHg by MerB, without the presence of Hg(II) reductive capability in the form of MerA, would be expected to be physiologically detrimental to cells due to an increase in reactive oxygen species as a result of a depletion of cellular redox buffers and increased enzyme destabilization by Hg(II) (Miller et al., 1991; Ariza and Williams, 1999; Valko et al., 2005). The observation of numerous genomes that encode MerB but not MerA and the observation that MerB in absence of MerR results in hypersensitivity to organomercury (Kiyono et al., 1995) suggests that cells must have evolved alternative mechanisms of mitigating the toxic effects of Hg(II) produced from MeHg demethylation. Here, potential alternative mechanisms of Hg(II) reduction in cells are examined with the goal of identifying putative functionalities that could mitigate the toxicity of Hg(II) generated by MerB in strains that lack MerA (Table 3). The suggested mechanisms could be tested for their role in organomercury detoxification to complement existing paradigms on microbial transformations in mercury biogeochemistry.

FIGURE 2

Venn diagrams reporting the number of genomes with MerA and/or MerB protein homologs. Analyzed genomes included those from Bacteria, Archaea, and plasmids that encode only MerA homologs, only MerB homologs, or that encode both MerAB. The number of genomes in each category is indicated.

TABLE 3

Possible Hg(II) detoxification mechanisms independent of the mer system1.

Several non-mer Hg(II) reduction mechanisms have been described [for a comprehensive list of such mechanisms see Grégoire and Poulain (2018)] and these, if coupled to organomercury degradation by MerB, could result in detoxification (Table 3). Constitutive Hg(II) reduction was observed in the dissimilatory metal-reducing bacteria Shewanella oneidensis MR-1, Geobacter sulfurreducens PCA, and Geobacter metallireducens GS-15 (Wiatrowski et al., 2006). Strain MR-1 does not encode MerA and while Geobacter do encode divergent MerA homologs, their direct role in Hg(II) reduction is not clear (Hu et al., 2013). For S. oneidensis MR-1 cells but not G. sulfurreducens PCA or G. metallireducens GS-15 cells, Hg(II) reduction was detected when they were provided with oxygen or fumarate as electron acceptors. However, all three strains were able to reduce Hg(II) when ferric oxyhydroxide was provided as a terminal electron acceptor (Wiatrowski et al., 2006). Subsequently, it was shown that the product of ferric oxyhydroxide reduction, magnetite, is likely responsible for the reduction of Hg(II) in these cells (Wiatrowski et al., 2009). This conclusion was further supported when Hg(II) was shown to be reduced by biogenic magnetite in the membrane of the magnetosome of two magnetotactic bacteria (Magnetospirillum gryphiswaldense MSR-1 and Magnetospirillum magnetotacticum MS-1) (Liu and Wiatrowski, 2018). Thus, Hg(II) reduction in these cells likely occurs via coupled iron and Hg(II) redox reactions as an alternative mechanism of Hg(II) reduction. This mechanism is therefore a slight variation on the reported reduction of Hg(II) by chemotrophic and mercury-resistant acidophilic ferrous iron-oxidizing thiobacilli (Sugio et al., 2003). In these cells, Hg(II), instead of oxygen, is reduced at the terminus of an electron transport chain that is initiated by the oxidation of ferrous iron.

Hg(II) can also be reduced by cellular processes that result in the production of excess reducing power (Table 3). For example, Hg(II) was reduced in cultures of purple non-sulfur bacteria (Rhodobacter capsulatus, Rhodobacter sphaeroides, and Rhodopseudomonas palustris) grown photoheterotrophically (Gregoire and Poulain, 2016). Using light as an energy source and acetate or butyrate as carbon source, these purple non-sulfur bacteria grow better in presence of sublethal concentrations of Hg(II) than without it. The authors elegantly show that by serving as an additional electron sink, Hg(II) can restore redox homeostasis thereby enhancing the overall metabolism of the cells (Gregoire and Poulain, 2016). Similarly, in cultures of obligate anaerobes from the order Clostridiales (Clostridium acetobutylicum ATCC 824 and Heliobacterium modesticaldum Ice1), MerA-independent Hg(II) reduction was observed during fermentation of pyruvate, presumably due to the production of reduced redox cofactors such as low potential ferredoxin (Grégoire et al., 2018). Similarly, H. modesticaldum Ice1 was shown to reduce Hg(II) when grown photoheterotrophically with pyruvate as a carbon source. In these cases, it is possible that low potential ferredoxins produced during phototrophic or organotrophic metabolism are responsible for Hg(II) reduction (Grégoire et al., 2018). In the examples outlined above, Hg(II) is reduced due to its electrophilic properties, leading to the production of volatile Hg(0).

Previous bioinformatics studies have shown that oxygenic Cyanobacteria do not encode MerA (Barkay et al., 2010), and the expanded MerA database reported here supports this conclusion with the exception of a single MerA-encoding Halothece genome (IMG ID: 2795385387) recovered from a halite pinnacle. Despite not encoding MerA, the cyanobacterium Synechocystis PCC6803 was shown to reduce Hg(II), likely mediated by a redox homeostasis system powered by low potential ferredoxins generated by photosystem II, like those mentioned above (Marteyn et al., 2009). Subsequently, a MerA-like enzyme was identified that was dependent on interactions with glutaredoxin for its activity; that activity resulted in the reduction and increased tolerance to Hg(II) as well as to uranyl acetate (Marteyn et al., 2013). While such mechanisms may help explain light-dependent Hg(II) reduction by phototrophs (Ben-Bassat and Mayer, 1978; Mason et al., 1995), none of the genomes that encode MerB and that lack homologs of MerA were derived from Cyanobacteria (Supplementary Table 1).

Finally, an array of studies have documented Hg(II) reduction by abiotic mechanisms that themselves are often indirectly connected to specific microbial activities and associated metabolites (Table 3). Dissolved organic matter (DOM) is an enormous effector of mercury bioavailability and this is largely attributable to interactions with organo-thiols (RSHg+) and Hg(II) (Ravichandran, 2004). Likewise, Hg(II) has been shown to be reduced by humic acids associated with natural organic matter (NOM) in aquatic environments (Gu et al., 2011). Reduction of Hg(II) by humic acids appears to be common in anoxic environments with abundant NOM explaining numerous observations of mer-independent formation of Hg(0) in such environments (see Grégoire et al., 2018). Interestingly, intricate interactions between microbes and NOM have been suggested to impact Hg(II) reduction. For example, NOM-enhanced Hg(II) reduction by S. oneidensis MR-1 was shown to occur by diverting electrons generated during anaerobic respiration to Hg(II) via NOM (Lee et al., 2018). Further, inhibition of Hg(II) reduction was observed at higher concentrations of NOM, likely due to sequestration of Hg(II) by RSHg+ groups in the NOM. Additionally, reduction of Hg(II) in the presence of DOC has been observed under sunlight in many different types of natural waters (Costa and Liss, 1999). This activity is affected by many factors such the chemistry and molecular weight of the DOC pool, the availability and wavelength of light, pH, dissolved oxygen, and the composition and abundance of cations and anions (Luo et al., 2020). Photoreduction of Hg(II) by intracellular algal components has also been documented (Kritee et al., 2018); however, the relevance of this process in mitigating the intracellular toxicity of Hg(II) produced from MeHg remains to be determined in species whose genomes encode MerB but not MerA.

Not surprisingly, several known alternative Hg resistance mechanisms are based on the interactions of Hg(II) with reduced sulfur (Table 3). Furthermore, uptake and transformations of Hg(II) in Hg(II)-methylating or reducing bacteria have been shown to be facilitated by thiols (-SH) and RSHg+ on proteins (Yu and Fein, 2017; Thomas et al., 2018, 2019; Song et al., 2020). Such interactions drive Hg behavior at all levels from the intracellular molecular level (Barkay et al., 2003) to the ecosystem level (Ravichandran, 2004) and they often lead to increased Hg(II) resistance. For example, the well characterized mer operon in the thermophilic bacterium Thermus thermophilus carries a gene, oah2, that codes for the synthesis of the low molecular weight thiol bacillithiol (Norambuena et al., 2018, 2019). Deletion mutants of oah2 as well as those of several genes involved in the reduction of bacillithiol were as sensitive to Hg(II) as a MerA deletion mutant. Furthermore, it was shown that bacillithiol can be oxidized by Hg(II). The authors suggested that along with Hg(II) reduction by MerA, bacillithiol buffers the cell against Hg(II) toxicity (Norambuena et al., 2018, 2019). In another example, the genes encoding capabilities that allow the obligate anaerobe Clostridium cochleariumm T-2P to demethylate MeHg and to generate hydrogen sulfide were encoded on a plasmid (Pan-Hou and Imura, 1981). Curing the strain of the plasmid resulted in a loss of both activities and this strain was more sensitive to mercurials than the wild type (Pan-Hou and Imura, 1981). This suggests that MeHg demethylation, likely mediated by MerB, is coupled to H2S generation and this results in the formation of insoluble mercury sulfide (HgS; cinnabar). A third example was described in two moderately thermophilic Bacillus and one Ureibacillus strains that were exposed to Hg(II) which resulted in the formation of a black precipitate, subsequently shown to be HgS by energy dispersive X-ray diffraction analysis (Glendinning et al., 2005). Volatile sulfide was not detected in the growth medium of these cells and it was therefore suggested that enhanced Hg tolerance in these cells is based on the production of non-volatile intracellular thiols that promote sequestration of Hg(II) as HgS (Glendinning et al., 2005). A better understanding of the role of thiols in Hg(II) sequestration in taxa whose genomes encode MerB but not MerA and targeted experiments where strains defective in thiol production are used as controls are needed to determine their potential role in mitigating the toxicity of Hg(II) generated by MeHg demethylation.

Many studies have documented additional interactions of Hg(II) with sulfur (Ravichandran, 2004), that might be considered for their role in Hg(II) resistance. For example, HgS complexes have been documented in Klebsiella pneumoniae M426 (Essa et al., 2002). This strain was shown to produce volatile organosulfur products that led to the formation of amorphous precipitates that contained Hg, carbon, and sulfur. Likewise, in some Cyanobacteria, intracellular deposits of HgS were observed and these were suggested to be formed via the precipitation of Hg(II) by intracellular thiol groups (Lefebvre et al., 2007).

The above mechanisms provide a framework for rationalizing how organisms that encode MerB, but not MerA, could demethylate Hg and mitigate the toxic consequences of intracellularly produced Hg(II). The mechanisms described are found in diverse taxa operating a variety of energy metabolisms ranging from aerobic to anaerobic and employing fermentative, heterotrophic, and phototrophic/photosynthetic metabolisms, indicating that MerA-independent Hg(II) detoxification can take place in a variety of ecological contexts. These mechanisms are possibly widespread in biology considering that thiols are universal components of all cells (Rooney, 2007; Benedikter et al., 2018), sulfide is produced in many cells through dissimilatory metabolism [e.g., sulfate or sulfur reduction (Plugge et al., 2011) and assimilatory metabolism; e.g., biosynthesis of cysteine (Kitabatake et al., 2000; Mukai et al., 2017)], and iron metabolism is common place among anaerobic or facultatively anaerobic cells (Weber et al., 2006). Further, reduced cofactors such as low potential ferredoxins are ubiquitous in anaerobic cells, including a variety of phototrophs and chemoautotrophs, the latter of which generate this low potential reductant through flavin-based electron bifurcation or ion-coupled translocation mechanisms (Buckel and Thauer, 2018; Boyd et al., 2020). Such low potential ferredoxins can also be generated through heterotrophic metabolisms, including those that respire, ferment, or photo-assimilate organic carbon (Gregoire and Poulain, 2016; Poudel et al., 2018). The fact that Hg(II) detoxification can be accomplished through these means, when combined with data presented herein revealing the diversity of taxa that encode MerA, indicates that MerA is a superior Hg(II) reducing catalyst with a dedicated role to reduce and thus detoxify Hg(II). The majority of the MerA-independent mechanisms of detoxifying Hg(II), as described above (Table 3), are potentially important to species inhabiting anoxic habitats. Consequently, the distribution of these mechanisms largely among anaerobes is consistent with the emergence of MerA among aerobes, which has been suggested to have taken place only after O2 began to accumulate in the biosphere and when HgS minerals were oxidatively weathered thereby increasing the release of toxic Hg(II) (Barkay et al., 2010; Boyd and Barkay, 2012).

Phylogenetic Reconstruction of MerA Homologs

The discovery of previously unidentified and diverse MerA protein homologs motivated a phylogenetic analysis to further evaluate the evolutionary history of MerA-mediated Hg(II) reduction (Figure 3). In the presented phylogenetic tree, MerA lineages were collapsed at the phylum level or at the class level in the case of Proteobacteria. Lineages were not collapsed when a monophyletic clade could not be forced at the taxonomic ranks outlined above or when MerA belonged to multiple phyla within a lineage. When a lineage was collapsed, it was additionally assigned a numeric designation (e.g., K02 Crenarchaeota; Figure 3; Supplementary Table 4) and the number of MerA homologs in each lineage is depicted in brackets. Non-collapsed clades were not provided with a “K” designation. Rapid bootstraps values were also calculated and are presented in Supplementary Dataset 1. Plasmid encoded MerA protein homologs were excluded from the analysis due to their polyphyletic nature (Hülter et al., 2017).

FIGURE 3

Phylogenetic tree of MerA protein homologs (plasmid MerA not included). Archaeal and bacterial clades are denoted with blue and red color, respectively, while purple-colored clades denote mixed lineages. Colors covering labels of lineages denote the three primary groups of MerA: Group A (light blue), Group B (light orange) and Group C (light green). The outer boxes right to the labels show the presence of MerB (light red box) and MerB-like (gray boxes) homologs in the genomes from which the MerA homologs were derived. Numbers in the light red and gray boxes denote how many of the MerA-encoding genomes within each collapsed node also encoded MerB or MerB-like homologs, respectively. The outgroup consists of 6 bacterial and archaeal dihydrolipoamide dehydrogenase sequences. Tree scale is in amino acid substitutions per site. Yellow stars denote lineages with representatives which their MerA activity has been confirmed. Single lineages not grouped with other clades are depicted with their gene ID and phylum assignment.

Phylogenetic reconstruction of MerA proteins showed that lineages broadly correspond to the phylum taxonomic level of the genome from which they originate, similar to prior analyses of MerA phylogeny (Barkay et al., 2010; Boyd and Barkay, 2012). The MerA phylogeny consists of three large empirically defined groups of sequences: Group A that consists exclusively of MerA from archaeal lineages, Group B that consists primarily of MerA from archaeal lineages along with several recently described bacterial lineages (i.e., Kryptonia, Candidatus Giovannonibacteria, Candidatus Peregrinibacteria among others), δ-Proteobacteria, and Chloroflexi, and Group C that consists primarily of MerA from bacterial lineages as well as MerA from a Euryarchaeota lineage (Figure 3; Supplementary Table 4). Group B and Group C include clades with MerA sequences from confirmed Hg(II) reducers (Boyd and Barkay, 2012; Supplementary Table 1). These clades have been denoted with yellow stars (Figure 1).

The earliest branching lineage of MerA (Group A) comprises only archaeal sequences, suggesting an origin for MerA among Archaea. Specifically, the archaeal phyla that comprise Group A include MerA from Crenarchaeota (class Thermoprotei; lineage K03), Thaumarchaeota (order Aigarchaeota; K02), Geothermarchaeota (lineage K04), Thaumarchaeota (K05), Euryarchaeota (class Thermoplasmata; K06 and K09), and Candidatus Geoarchaeota (K07 and K10). Characterized members of these phyla are thermophiles and often acidophiles and have been identified in a variety of thermal settings including hydrothermal fluids, geothermal springs, and other volcanic-influenced environments (Golyshina, 2011; Kozubal et al., 2013; Reysenbach, 2015; Zillig and Reysenbach, 2015; Jungbluth et al., 2017; Hua et al., 2018). Of the above lineages, MerA activity has been documented in the crenarchaeota Sulfolobus solfataricus (Schelert et al., 2004) and Metallosphaera sedula (Artz et al., 2015). Apart from the lineages Crenarchaeota and Euryarchaeota that include homologs from aerobic archaea (Golyshina et al., 2009; Kozubal et al., 2011), all thermophilic and hyperthermophilic archaea with MerA that cluster in Group A are putative anaerobes. The most basal lineage comprises both aerobic or facultatively anaerobic crenarchaeota, including MerA from the genera Pyrobaculum, Thermoproteus, Thermofilum, and Aeropyrum. This suggests a thermophilic origin for MerA in a volcanically influenced habitat, potentially with low O2 availability, a finding that is consistent with previous metagenomic surveys (Geesey et al., 2016), phylogenetic studies of MerA (Barkay et al., 2010; Boyd and Barkay, 2012), and the unusually high optimal temperature for MerA enzymes in a variety of mesophilic bacterial lineages (Vetriani et al., 2005; Barkay et al., unpublished). Volcanic habitats are often enriched in naturally high levels of Hg(II) (Varekamp and Buseck, 1986), possibly due to oxidative leaching of Hg-containing minerals (e.g., HgS) in the subsurface. Naturally high levels of Hg(II) in geothermal environments would be expected to be a potent selective pressure to evolve a dedicated Hg(II) detoxification system in the form of MerA. Thus, the present phylogeny, which consists of a much larger data base of MerA sequences (Table 1) than in previously analyses (Barkay et al., 2010; Boyd and Barkay, 2012), has resulted in a refinement of the MerA evolutionary path to firmly establish an origin for MerA among thermophilic Archaea from geothermal environments.

MerA homologs that comprise Group B consist of both archaeal and bacterial homologs that are further separated into 3 subgroups denoted by light blue and purple coloring (Figure 3). The first subgroup consists of MerA homologs from thermophilic and hyperthermophilic euryarchaeal classes Thermococci and Archaeoglobi (K11). This MerA lineage is sister to a lineage that comprises MerA from members of Candidatus Daviesbacteria (K12), Candidatus Peregrinibacteria (K13), Candidatus Giovannonibacteria (K14), and a variety of bacterial MerA homologs that belong to members of the Candidate Phyla Radiation (CPR) recovered from an aquifer close to the Colorado River, United States (Brown et al., 2015; Anantharaman et al., 2016). This suggests genomic acquisition of MerA among a bacterial member of the CPR through lateral gene transfer (LGT) that likely involved a thermophilic archaeon. The second subgroup of MerA homologs in Group B includes lineages from clades K15 to K21 (Figure 3; Supplementary Table 4). The second subgroup of MerA homologs in Group B includes all lineages from clades K15 to K21 (Figure 3; Supplementary Table 4). Clade K15 consists of MerA homologs from euryarchaeal methanogenic and methane-oxidizing Methanomicrobia (Kendall and Boone, 2006; Arshad et al., 2015), while the Thaumarchaeota clade K16 includes homologs from aerobic ammonia-oxidizing Nitrososphaera (Kerou and Schleper, 2016). The remaining MerA homologs in this subgroup (K17-K21) are from thermophilic Thermoprotei and Thermoplasmata. The third subgroup in Group B (depicted by orange) is nested within the second archaeal subgroup and includes many bacterial MerA such as the putatively anaerobic gammaproteobacterial group K27 (Figure 3; Magot et al., 1992; Vandieken et al., 2006; Shelobolina et al., 2007; Slobodkina et al., 2012), the anaerobic Kryptonia, isolated from neutral pH thermal springs (Eloe-Fadrosh et al., 2016), the anaerobic Chloroflexi from the Dehalococcoides genus (K28; Löffler et al., 2013), and the anaerobic Atribacteria (Liu et al., 2019). Acquisition of MerA among Bacteria comprising the third subgroup of Group B likely occurred through a LGT event involving an archaeon, possibly an ancestor of Candidatus Odinarchaeota, Candidatus Bathyarcheota, or Candidatus Heimdallarchaeota. In summary, multiple cross-domain LGT events are invoked to explain the genomic distribution of MerA that comprise Group B, events that substantially increased the taxonomic and ecological distribution of MerA in natural systems.

Group C of the MerA phylogeny comprises homologs that are primarily bacterial, with the exception of one euryarchaeal lineage (K34; Figure 3). The first subgroup in Group C (colored red; Figure 3) includes MerA from the thermophilic bacterium Deferribacter desulfuricans (Takaki et al., 2010) and several MerA homologs from unclassified bacterial lineages (lineage K32) that are reported as microaerophilic thermophiles according to the IMG database. This lineage also includes several MerA sequences from MAGs that are derived from artificial or human-impacted environments (Supplementary Table 2). These include the manganese-oxidizing Nitrospiraceae bacterium Mn-1 isolated from tap water (Yu and Leadbetter, 2020), Tissierella creatinophila KRE4 isolated from sewage sludge (Harms et al., 1998), several unclassified deltaproteobacteria (lineage K31) isolated from acid mine drainage (Tan et al., 2019), and several unclassified deltaproteobacteria and unclassified Nitrospirae bacteria (lineages K29 and K30, respectively) from a metal-contaminated aquifer near the Colorado River, Colorado, United States (Hug et al., 2015). This subgroup of Group C also includes MerA from the thermophilic hydrogen- or sulfur-oxidizing lineages of Aquificae (Gupta, 2014), a bacterial lineage where MerA was previously suggested to have originated (Barkay et al., 2010; Boyd and Barkay, 2012). However, Aquificae MerA are nested among several bacterial lineages including unclassified deltaproteobacterial (K29 and K31) and an unclassified Nitrospira MerA lineage (K30), all of which derive from MAGs. The discrepancy between the findings of previous studies that showed Aquificae MerA to branch at the base of MerA phylogenies (Barkay et al., 2010; Boyd and Barkay, 2012) and the findings reported herein is almost certainly due to the nearly ten-fold increase in the diversity of archaeal MerA considered in the present study, including those from MAGs and SCGs from environmental samples.

A second subgroup within Group C (colored purple, Figure 3) include MerA from two Archaeoglobus (lineage K34; Euryarchaeota) genomes from deep basaltic subseafloor fluids collected from the Juan de Fuca ridge (Jungbluth et al., 2017). Within this subgroup, lineage K34 branches basal to MerA from several Chloroflexi lineages (lineages K35, K36, K39), a Firmicutes bacterium (Limnochorda pilosa HC45), Nitrospirae lineage (K37), a lineage comprising unclassified Bacteria (K38), and the Acidobacteria (lineage K40). Apart from several MerA-encoding Chloroflexi genomes that derive from anaerobic aquatic sediments (Löffler et al., 2013; Hug et al., 2015), the rest of the MerA-encoding Chloroflexi lineages are putative thermophiles from geothermal environments (Hedlund et al., 2015; Watanabe et al., 2016; Sakoula et al., 2018; Zheng et al., 2019). The position of Euryarchaeota MerA in the Group C lineage suggests at least two separate cross-domain LGT events possibly took place to genomically acquire MerA as well as to diversify it among members of Group C, with both likely taking place in geothermal environments. The first LGT event likely involved genomic acquisition of MerA in a thermophilic ancestor of Euryarchaeota involving what was likely a thermophilic bacterium from the first subgroup (purple) of Group C, such as an ancestor of the Deferribacteres. The second LGT likely involved transfer of MerA from a thermophilic ancestor of Euryarchaeota to a thermophilic ancestor of Chloroflexi.

The third subgroup of Group C MerA in the phylogeny (red, Figure 3) comprises two predominant lineages of bacteria, both of which include early evolving members of the Proteobacteria (within this subgroup of Group C). Thus, it is possible that MerA in this subgroup was genomically acquired via LGT involving an ancestor of Proteobacteria. Interspersed among proteobacterial MerA are MerA from a variety of bacterial phyla including lineages comprising non-thermophilic members of the Proteobacteria, Actinobacteria, Bacteroidetes, Verrucomicrobia, Chloroflexi, Candidatus Rokubacteria, among others. In addition, this subgroup of Group C includes MerA from both thermophilic and non-thermophilic members of the Deinococcus-Thermus lineage. This suggests numerous LGT events have taken place among ancestors of these phyla that have broadened both the taxonomic diversity and ecological distribution of MerA. This is consistent with the existence of selective pressure to acquire and maintain MerA in taxa that inhabit a variety of environment types and supports previous studies that have shown that mer genes have a high propensity to be laterally transferred (Bogdanova et al., 1998; Barkay et al., 2003; Ojo et al., 2004), including several instances of cross-domain transfers (Barkay et al., 2010; Boyd and Barkay, 2012).

Mapping the Taxonomic Distribution of MerB Onto the MerA Phylogeny

Although proteins with TRASH domains (a metal-binding domain involved in heavy metal sensing, trafficking, and resistance) and the NosL protein (accessory protein of the nitrous oxide reductase system) have been suggested to share ancestry with MerB due to structural homology (Ettema et al., 2003; Taubner et al., 2006; Kaur and Subramanian, 2014), alignments revealed very little overlap in or homology among aligned positions (data not shown). As such, these paralogs cannot confidently be used as an outgroup in phylogenetic reconstructions of MerB (Supplementary Figure 2). Other suitable paralogs do not exist in current databases, preventing a directional assessment of the evolution of this protein. As an alternative, the distribution of MerB among MerA-encoding genomes was mapped on the MerA phylogeny (Figure 3). Assuming that MerA is the backbone of the mer detoxification system and that other functionalities were recruited later (Boyd and Barkay, 2012), this approach can provide insight into when MerB was recruited to MerA-encoding genomes as a dedicated strategy to detoxify organomercurials and increase the spectrum of mercurial resistance. Importantly, concentrations of MeHg and other organomercury compounds are often below the toxicity limit of cells in most environments. For example, MeHg concentrations even in contaminated environments are at ppt to ppb range (Gilmour et al., 1998) while MeHg toxicity to microbial growth is at the ppb to ppm range (Baldi et al., 1993). This raises the question of whether such compounds served as agents of natural selection as it related to acquisition of MerB and, if not, what the selective pressure to acquire MerB was.

MerB proteins were absent among members of the early branching lineages of crenarchaeal and thaumarchaeal MerA (Group A, Figure 3). Additionally, MerB homologs were not identified among members of the Candidatus Geoarchaeota, an archaeal phylum that features multiple genomes that encode MerA homologs (Figure 3; Table 1). The earliest branching MerA homolog in a genome that also encoded MerB was identified in Group B and specifically in a subgroup that includes lineages K15 and K16 (Figure 3; Supplementary Table 4). Lineage K15 includes methanogenic and mesophilic methane-oxidizing Euryarchaeota (Arshad et al., 2015), and lineage K16 includes mesophilic ammonia-oxidizing Thaumarchaeota (Kerou and Schleper, 2016). This suggests that MerB was likely recruited to the mer operon in an ancestor of these phyla. Previous findings revealed an increase in mer-encoded functionalities over evolutionary time (Boyd and Barkay, 2012). Similarly, early evolving lineages such as Crenarchaeota, Thaumarchaeota, and Geothermarchaeota tend to feature only MerA, while more recently evolved lineages (e.g., Proteobacteria) tend to harbor both MerA and MerB and other ancillary functions to improve Hg tolerance. This may suggest that recruitment of MerB in MerA-encoding organisms occurred after they had diversified into more moderate temperature environments.

Other Group B MerA-encoding lineages that also encode MerB include a member of the archaeal phylum Heimdallarchaeota, a phylum of organisms that have been detected in anaerobic marine sediments (Zaremba-Niedzwiedzka et al., 2017). Additionally in Group C, MerA-encoding genomes also often encode MerB, including a MerA-encoding Nitrospiraceae bacterium isolated from municipal tap water and a MerA-encoding Firmicutes bacterium isolated from sewage sludge (Harms et al., 1998). Many of the natural and artificial environments where organisms encoding both MerA and MerB have been isolated from tend to also harbor microorganisms that encode hgcAB genes (Podar et al., 2015; McDaniel et al., 2020). The production of MeHg in these habitats by HgcAB may have been the selective pressure that led to the recruitment of MerB to these organisms allowing for MeHg detoxification. The low (nM) concentrations of MeHg and other organomercurials in most natural habitats (Gilmour et al., 1998; Schartup et al., 2015) may reflect a steady state concentration attributable to near balanced rates of production and degradation.

Early evolving subgroups of Group C MerA include thermophilic Chloroflexi lineages; specifically, members of the Ktedonobacteria class of Chloroflexi (lineage K36) from geothermal hot springs (Tomazini et al., 2019) and Kallotenuaceae members (lineage K39) from the Great Boiling Spring in Nevada (Hedlund et al., 2015). However, these MerA-encoding genomes tend to not encode MerB. High temperature environments like hydrothermal vents and geothermal springs apparently do not favor microbial Hg(II) methylation, based on the limited distribution of hgcAB in these environment types (Driscoll et al., 2013; Podar et al., 2015). Despite reports of hgcAB, and thus Hg(II) methylation potential, in metagenomes recovered from geothermal environments, none of the organisms that encode these HgcAB can be confidently assigned to thermophiles. For example, while hgcAB was identified in marine mats in the submarine Kolumbo volcano (Podar et al., 2015), these mats were collected from a low-temperature zone in this system (Oulas et al., 2015). Furthermore, fused hgcAB sequences were identified in the thermophilic archaeon Pyrococcus furiosus (Podar et al., 2015). However, this strain was shown to be incapable of methylating Hg(II). Additionally, BLASTp searches of hgcA identified a homolog of the gene in a member of the archaeal Thermoplasmatales. However, this microorganism is an anaerobic mesophile isolated from an aquatic environment (Baker et al., 2015). Finally, hgcA was identified in the Coriobacteriia class of Actinobacteria and was linked to a geothermal environment because the Coriobacteriia belongs to the OPB41 class (McDaniel et al., 2020). The OPB41 class was first characterized through 16S rRNA gene sequencing of DNA from a hot spring (Hugenholtz et al., 1998). However, members of OPB41 have been identified in numerous non-thermophilic environments, including in cold basaltic sea sediments (Bird et al., 2019). As such, it is unclear if this hgcAB-encoding organism is a thermophile.

Collectively, these observations suggest that the aforementioned lineages that encode both MerA and MerB likely recruited MerB in a lower temperature habitat where mercury methylation is more likely to occur (Celo et al., 2006; Martins et al., 2006). Organisms with MerB then radiated among numerous bacterial and archaeal phyla that inhabit natural and anthropogenically impacted environments where Hg(II) is available for anaerobic methylation. On the basis of the distribution and abundance of hgcAB genes, such environments predominantly include subsurface aquifer systems, thawing permafrost layers, aquatic sediments, wetlands, and saturated soils (Podar et al., 2015; McDaniel et al., 2020).

Conclusion

The present study, which included analyses of the distribution and diversity of MerAB in isolate genomes, SCGs, and MAGs, significantly expanded the taxonomic and phylogenetic diversity of homologs when compared to previous studies that analyzed only the genomes of isolates (Barkay et al., 2010; Boyd and Barkay, 2012). This expansion included identifying MerA and/or MerB in genomes that correspond to 55 bacterial and archaeal lineages at the phylum level, including numerous phyla that were not previously reported to encode these functionalities. The expansion of the diversity of taxa that carry the mer system clearly points to the need to characterize the Hg transformation activities of a larger and more diverse range of organisms. The insights obtained will improve our understanding of the physiology and ecology of microorganisms that largely modulate the forms and fate of Hg in the environment, their evolution, and their potential uses in Hg bioremediation.

Numerous microbial genomes were identified that encode MerB but not MerA, suggesting the presence of potentially alternative mechanisms to detoxify cytoplasmic Hg(II) formed from demethylation of MeHg. A number of potential mechanisms that allow for Hg(II) detoxification are presented. These are largely based on alternative mechanisms for reduction of Hg(II) and on various aspects of the interactions of Hg with sulfur. Furthermore, the existence of these alternative mechanisms of detoxifying Hg(II) suggest that MerA originated as a superior Hg(II)-reducing catalyst that has been maintained in a diverse array of microbial taxa. The numerous successful LGT events required to explain the widespread taxonomic distribution of MerA across Archaea and Bacteria support this hypothesis. The broad distribution and high occurrence of MerB in the absence of MerA highlights how little we know about this enzyme and its role as a catalyst in the mer system. Moreover, in light of the scarcity of MerB paralogs, untangling the path of its evolution is challenging. Further, it is not clear what selective pressure might have shaped the origin and evolution of MerB considering the low concentrations of MeHg in the natural environment (although, that concentration might be higher immediately proximal to methylating microbes – see below). Thus, further exploration of MerB homologs, including those with sequence variations of unknown functional consequences, in diverse taxa from a variety of environmental settings (e.g., low versus high temperature) is warranted.

With an expanded database of MerA homologs and a robust phylogenetic reconstruction of these homologs, we present a refined path for the evolution of MerA that improves our understanding of the origin, evolution, and environmental distribution of this enzyme. Firstly, the earliest evolving MerA proteins were identified among the genomes of organisms inhabiting high temperature hydrothermal environments, which is consistent with previous studies (Barkay et al., 2010; Boyd and Barkay, 2012). It is possible that naturally elevated concentrations of Hg(II) that typify hydrothermal environments acted as the original selective pressure to evolve this enzyme. However, the earliest branching MerA sequences were from thermophilic Archaea rather than Bacteria, which differs from the origin for MerA among bacterial and thermophilic Aquificales as advocated previously (Barkay et al., 2010; Boyd and Barkay, 2012). The cultivation-independent nature of the SCGs and MAGs included in the study reported here, which resulted in a 25-fold increase in the number of MerA homologs known, is almost certainly responsible for the change in the topology of the tree. Nonetheless, like previous phylogenetic studies of MerA, numerous instances of cross domain and cross phylum LGT events were observed in the phylogenetic reconstruction conducted herein. Combined with the finding of multiple (up to 8 copies) of MerA in the genomes of some organisms and the detection of MerA on plasmids (which are likely to be underreported for reasons outlined above), these data point to the fitness advantage that MerA-encoding cells have in Hg(II) contaminated environments and the selective pressure for cells to maintain MerA once it is acquired via LGT. The duplication of MerA encoding genes may have allowed for cells to increase their responsiveness to Hg(II) toxicity (gene dosage effect) or to respond to Hg(II) toxicity with an enzyme with slightly altered kinetics or substrate specificity (sub-functionalization).

While it was not possible to perform an analysis of the evolution of MerB using a rooted phylogeny, mapping of MerB on the MerA phylogeny suggests that it was likely recruited to the mer operon relatively recently. Furthermore, the near universal absence of MerB homologs among thermophiles, when combined with the fact that the earliest evolving MerA-encoding genome that also encodes MerB derives from a mesophile, suggests that the ability to demethylate MeHg emerged in a lower temperature setting following the emergence of mer systems in mesophiles. This is potentially consistent with the limited distribution of genes involved in biotic MeHg production (i.e., hgcAB) in thermophiles or high temperature habitats but being pervasive in lower temperature habitats (Podar et al., 2015; McDaniel et al., 2020). Consistent with this notion, MerB-encoding genomes were identified among organisms or DNA recovered from a variety of habitat types that likely favor MeHg production, including aquatic sediments, thawing permafrosts, subsurface aquifers, and soils. This may suggest a relationship between the biological production and consumption of MeHg in natural habitats. The expanded taxonomic diversity of MerB homologs and the metabolic backgrounds of the organisms that they derive from may provide new strategies that could be employed in the management of MeHg contamination and its bioremediation.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author Contributions

CAC and ESB designed the study. CAC carried out bioinformatic and phylogenetic analyses. All authors contributed to data analysis and the writing of the manuscript.

Conflict of Interest

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.