Microbial Mercury Methylation in Winooski River and Englesby Brook, Vermont

1. Introduction

Mercury (Hg) is a toxic environmental contaminant that bioaccumulates in the food web, posing risks to both ecosystems and human health [1]. Major anthropogenic sources of Hg include industrial activities, coal burning, and gold mining [1,2]. Over the past several decades, studies have highlighted mercury’s detrimental effects on neurodevelopment, cognition, and cardiovascular health [3]. Furthermore, Hg interacts synergistically with other contaminants such as polychlorinated biphenyls (PCBs) and other metals, amplifying its toxic impact on organisms [2,4].

Methylmercury (MeHg), produced via mercury methylation, is the only form of Hg that biomagnifies in the food web. Unlike Hg(0), Hg(II), and Me₂Hg, which are non-reactive and easily diffuse, MeHg is more readily incorporated into biota [1]. Mercury methylation occurs in aquatic ecosystems, including wetlands, sediments, benthic surfaces, pelagic environments, and flooded or inundated areas like rice paddies [5,6,7]. This process is most prominent at oxic-anoxic boundaries [8,9,10].

In microbes, mercury methylation serves as a detoxification mechanism, enabling microbes to resist Hg toxicity [11,12,13]. Methylating microbes, previously identified as mainly facultative and obligate anaerobes such as sulfate-reducing and iron-reducing bacteria, are responsible for biotic mercury methylation [7,14]. Microbial uptake and release of MeHg is rapid, unlike higher organisms, which tend to bioaccumulate it [15]. Sulfate-reducing bacteria, in particular, have been identified as the primary methylators in environments such as estuary sediments [16]. These microbes likely utilize a corrinoid-containing sulfur or iron protein that acts as a methyl group donor in the acetyl-CoA pathway in the cytoplasm, and thus methylmercury synthesis likely competes with acetyl-CoA for methyl groups [17].

The hgcA and hgcB genes have been identified as essential for mercury methylation in Desulfovibrio desulfuricans ND132 and Geobacter sulfurreducens PCA [14]. The proteins encoded by these genes facilitate the methylation process by producing the corrinoid protein carrying out methyl transfer (HgcA) and creating 2[4Fe-4S] ferredoxin for corrinoid reduction (HgcB). Deletion of these genes reduced methylmercury production by over 99%, and reintroduction restored methylation activity. These genes are predominantly found in the phyla Proteobacteria and Firmicutes and across other bacterial and archael strains, such as Deltaproteobacteria, Clostridia, and Methanomicrobia [14,18]. In Geobacter sulfurreducens PCA, deletion of the hgcAB genes caused minimal growth disruption but altered metal ion homeostasis [19]. The presence of the hgcAB genes can predict Hg-methylating potential, as demonstrated by their detection in 15 previously untested bacteria with hgcAB homologues [20].

Uncertainties persist regarding the effects of climate change and land use on mercury cycling and methylation. Climate change, particularly temperature increases and altered hydrological fluxes due to flooding and drought, may exacerbate mercury methylation and bioaccumulation [2,7]. Predicted to intensify due to climate change, flooding and drought can trigger redox oscillations, enhancing methylation by sulfate-reducing bacteria through increased sulfide oxidation [21,22,23]. In the Cottage Grove Reservoir, Oregon, USA seasonally inundated sediments had significantly higher mercury methylation activity compared to permanently inundated areas [21]. Another study of the Upper St. Lawrence River, New York, USA wetland sediments estimated significant Hg mobilization from a government plan to restore river fluctuation to pre-Moses-Saunders hydropower dam oscillations, due to the presence of mercury-methylating microbes [5]. Intense flooding in Vermont in July of 2023 may play a role in THg load in sediments sampled for this study [24].

In the Lake Champlain Basin of Vermont, USA, Malletts Bay has been identified as an area with high concentrations of total mercury (THg) and MeHg [25]. The Winooski River and Englesby Brook have also been notable for relatively high Hg concentrations and were classified as impaired under Vermont’s Water Quality Standards for mercury [26,27,28]. Both sites were heavily flooded in July 2023, which likely intensified mercury mobilization and methylation in these areas, marking them as important sites of investigation [29]. This study seeks to uncover the potential consequences of legacy mercury deposition on microbial mercury methylation in Winooski River and Englesby Brook. Sediment samples from different zones (Riparian Land, RL; Riparian Water, RW; Bay, B) were analyzed to understand community composition and the likelihood of mercury methylation. To identify Hg-methylating microbes and quantify hgcA genes, we employed 16S rRNA amplicon sequencing and metagenomic analysis, respectively, as described by previous studies [30,31,32] and followed previously published sediment sampling and elemental analysis protocols [5].

2. Materials and Methods

2.1. Experimental Design

Samples were collected at two primary sites in Vermont: Winooski River in Colchester on October 28, 2023 (44.524444, -73.256389) and Englesby Brook in Burlington on November 4, 2023 (44.457778, -73.219722). At both Winooski River and Englesby Brook, there were three separate locations: the sediment within the inundated river (Riparian Water, RW), the riparian land at the edge of the river not covered at the time of sampling by water (Riparian Land, RL), and the bay area (Winooski Bay, WB and Englesby Bay, EB), which is the river mouth that feeds into Lake Champlain (Figure 1). For WB we sampled at Delta Park (44.5375323069514, -73.27796273800251). Meanwhile, we sampled for EB at Blanchard Beach (44.45819838869913, -73.2240284237561). From the RW and RL regions of each site, three cores from three plots of 1 m² were collected (Downstream 1, D1; Downstream 2, D2; Upstream, U) (Figure 1), whereas three cores in total were collected from a 1 m² plot at each bay area.

2.2. Sample Collection

Sediment and hydric soil cores were taken with a polycarbonate tube (1.5 m long, 2.54 cm external diameter, 2 cm internal diameter), hammering with a rubber hand mallet when necessary. After reaching a depth of around 10 cm, a vacuum was created inside the tubing by placing one hand tightly on top of the tube, then pulling the corer straight out of the ground and immediately transferring to a clean and labeled air-tight polyethylene bag by releasing the sealing hand. If the sediment was particularly loose, a 3D-printed plug was inserted at the bottom of the tube to avoid loss of sediment. Prior to sample collection, corers were acid-washed with 3% HCl for 15 minutes and DI water for 15 minutes. At the field sites, soil corers were equilibrated by rinsing with river water three times prior to taking the samples as well as between the samples. Samples were stored in a cooler during transport to the lab until further processing.

2.3. Sample Processing and Storage

Excess water was removed from polyethylene bags after the samples had settled and a clear separation between the sediment and water was visible. Within 1-3 days of sampling, samples were cleared of large debris (>3 mm) by sifting through using stainless-steel tweezers and homogenized using a stainless-steel coffee grinder (Cuisinart Grind Central; Stamford, CT) or acid-washed mortar and pestle, transferred to polypropylene tubes, then stored frozen at -20° C until further analysis. Subsamples of 1 g each were removed from homogenized sediment samples for carbon and nitrogen analysis. From each site in each location (D1 RL, D1 RW, D2 RL, D2 RW, U RL, U RW, and B for both Winooski (W) and Englesby (E)) subsamples of 2 g from each of the three cores were removed from thawed homogenized sediment samples, then pooled and homogenized again for DNA extraction.

2.4. DNA Extraction and Sequencing

DNA was extracted using the Zymo Research Soil Microbe DNA Microprep kit, following the protocols by the manufacturer (Zymo Research, Irvine, CA, USA). Duplicate extractions were performed for each of the 15 samples. Briefly, 0.25 g of each homogenized sample was lysed in a bead-beating tube, centrifuged, filtered through DNA-binding spin columns, eluted with nuclease-free water, and then transferred to a clean 1.5 mL microcentrifuge tube after filtering through a Spin II-microHRC Filter. A negative control of 0.25 mL of sterile nuclease-free H₂O was also extracted through the kit [35]. A Quant-iT™ PicoGreen™ dsDNA Assay (Thermo Fisher Scientific,Waltham, MA, USA) was conducted following DNA extraction to quantify the amount of DNA in each sample and the negative control. All 14 extracted DNA samples and one negative control were sent to Novogene (Sacramento, CA, USA) for Illumina 16S rRNA amplicon library prep and Illumina Miseq PE250 sequencing. Sequenced data were sent back as raw reads, then processed and analyzed using the DADA2 pipeline. Three extracted DNA samples (D1, D2, U) from each site (Englesby RL, Englesby RW, Englesby Bay, Winooski RL, Winooski RW, Winooski B) were pooled and shipped to Novogene for NovoSeq PE150 metagenomic sequencing.

2.5. Mercury Analysis

Total mercury content in the sediment samples was quantified in approximately 1 g of wet sediment from each core, using a Milestone (Sorisole, Italy) Direct Mercury Analyzer (DMA-80 evo). The instrument was initially calibrated using aqueous mercury standards and verified with duplicate measurements of IAEA-456: Coastal Sediment (International Atomic Energy Agency). Samples were dried at 200 C for one and a half minutes and combusted at 650°C for three minutes; the instrument was purged for 75 seconds and the amalgamator heated for 12 seconds within a 30 second signal recording window. Samples were analyzed in duplicate with relative standard deviations of less than 10%. In addition to IAEA-456, an instrument blank, ultrapure water blank, a procedural blank, and a matrix-spike were analyzed every 10 samples (20 boats). Acceptable percent recovery for a matrix spike was between 80-120%. The limit of detection was defined as three standard deviations above the mean of ultrapure water blanks, and the limit of quantification as defined as ten standard deviations above the mean. All samples were above the limit of quantification.

2.6. Carbon & Nitrogen (CN) Analysis

Subsamples (~1 g) of each sediment sample were freeze-dried and then ground with a mortar and pestle. A standard curve was created using four masses of aspartic acid: ~0.5 mg, ~2.5 mg, ~4.5 mg, and ~6.5 mg. These resulted in masses of nitrogen of 0.041 mg, 0.253 mg, 0.458 mg, and 0.672 mg, and masses of carbon of 0.171 mg, 0.917 mg, 1.617 mg, and 2.337 mg. The area of the lowest calibration for nitrogen was 90,234.5 mV*sec, while the area of lowest calibration for carbon was 872,005 mV*sec. A standard reference material, sCNS, was used during each analysis to confirm that the standard curve was accurate. Around 30.0 mg of each sample were analyzed using a Flash 2000 CNS-analyzer (Thermo Fisher Scientific, Waltham, MA, USA).

2.7. Data Analysis of Elemental Concentrations

All statistical analyses were conducted using Microsoft Excel and R Studio software (v. 4.3.1). The elemental concentration value from the three cores of each 1 m² plot (Figure 1) were averaged prior to analysis using dplyr (v. 1.1.4) [36]. Data normality and homogeneity was assessed using the Shapiro Wilk test stats package (v. 4.3.3) in base R [37]. The car package (v. 3.1.2) was used to perform a Levene test [38]. An ANOVA test was used to compare elemental concentrations between the sites, followed by a post-hoc Tukey’s (stats package, v. 4.3.3). The correlation between different elements (carbon and nitrogen) and the concentration of total mercury at each site was determined using a Pearson’s product-moment correlation test using the stats package (v. 4.3.3) [37].

2.8. 16S rRNA Amplicon Analysis

DADA2 (v. 1.6) was used to demultiplex, quality control (QC), merge paired-end reads, and assign amplicon sequence variants (ASVs) and taxonomic assignments to all samples and the blank [39]. After creating quality profiles for each forward and reverse .fastq sequence, we trimmed from 150 nucleotides to the end for the forward reads, and from 220 nucleotides to the end for the reverse reads. An ASV table was created using dereplicated unique sequences from merged forward and reverse reads. Chimeras were removed, and taxonomy was assigned to the species level using DECIPHER (v. 2.20.0), referencing the SILVA_SSU_r138 training set [40]. Taxonomic composition and diversity of samples were visualized using phyloseq (v. 1.42.0)[41]. Non-metric dimensional scaling (NMDS) using ASV relative abundance Bray-Curtis dissimilarity distance was used to compare dissimilarity in microbial community composition among the different sampling locations (RL, RW, Bay) in Winooski River and Englesby Brook. Differences between the Hg-methylating microbial community composition between the different sites were compared by an ANOVA test followed by a post-hoc Tukey’s test. Previously identified Hg-methylating taxa were filtered for and analyzed further [42]. Differences between the Hg-methylating microbial community composition between the different sites were compared by an ANOVA test followed by a post-hoc Tukey’s test.

2.9. Metagenome Analysis

Raw metagenomic reads were quality filtered and trimmed using trimmomatic (v 0.39) [43], then assembled using megahit (v1.2.9) [44]. Open reading frames (ORFs) within the assembled contigs were identified using prodigal (v2.6.3) [45] and converted to an amino acid sequence file (.faa). HMMER (v3.4) [46] was used to search for sequences that matched the Oak Ridge National Lab (ORNL) Hg-methylating database [47]. Following the ORNL README file, we placed aligned sequences onto a maximum likelihood reference HgcA tree and assigned taxonomy to sequences based on their placement on the tree. Contigs identified via HMMER as HgcA sequences were also classified using the Basic Local Alignment Search Tool (BLAST) [48], and the top hits were recorded.

3. Results

3.1. Carbon, Nitrogen, and Hg Analyses

Although average THg concentrations (ng/g) were low overall (~10-40 ng/g), the Englesby Brook Riparian Land (ERL) region had greater average total Hg (THg) concentration (ng/g) than any of the other sampling sites in both Winooski River and Englesby Brook (Figure 2). Englesby Brook samples also had greater average THg concentrations (41.93 ng/g) overall compared to Winooski River. In Winooski River, average THg concentrations were significantly positively correlated with average carbon concentrations (mg/g) (p=0.02415), but not with average nitrogen concentrations (mg/g) (Table A1, Figure 3). However, there was no significant relationship between average THg and either average carbon or nitrogen in Englesby Brook (Table A2).

There was no significant difference in carbon or nitrogen mean concentrations between Winooski River and Englesby Brook (respectively, p=0.1009, p=0.2223) (Figure A2.1; Figure A2.2). Winooski Bay had a mean carbon concentration of 36 mg/g dry weight, significantly greater than both Winooski RL (4.5 mg/g dry weight) and Winooski RW (6.9 mg/g dry weight) (p=0.0264) (Table A1, Figure A2.1A). However, there was no significant difference in mean nitrogen concentrations between sampling sites in Winooski River (p=0.3733) (Figure A2.2A). Meanwhile, for Englesby Brook there was no significant difference in the mean carbon (p=0.4701) or nitrogen concentrations (p=0.4273) between the sampling sites (Table A2, Figure A2.1B, Figure A2.1B). Overall, nitrogen concentrations were low, and many were below detection (detection limit = 0 mg/g, calculated as [mean of blanks + standard deviation +/-3]) for all sampling sites in both Winooski River and Englesby Brook).

3.2. Microbial Community Structure Based on 16S rRNA

Non-metric multidimensional scaling (NMDS) analysis of 16S rRNA data showed microbial community differences between Winooski River and Englesby Brook (Figure 4). No clear trends by sampling site were observed, though Englesby Brook Riparian Land samples clustered together, while others were scattered. Notably, the Winooski Riparian Water sample was distant from other Winooski samples, likely due to the dominance of Gammaproteobacteria abundance compared to other Winooski samples (Figure 5).

A high diversity of microbes was identified at the phylum level at all sites. Proteobacteria was the most abundant phylum across all 10 sites (Figure A2.3). Other abundant phyla included Actinobacteria, Firmicutes, and Bacteroidota. One Winooski RW site, WD1 RW, had the least observed phylum-level diversity for the 16S rRNA data sets. At the class level, the most abundant taxa across all 10 sites were Gammaproteobacteria, Desulfomonilia, and Bacilli (Figure 5). Overall, genera identified as being capable of Hg-methylation were potentially quite high in abundance in all sampling sites, ranging from 0-5% of ASVs within the samples (Figure 6). Winooski D1 RL and Englesby URL had the greatest abundance of taxa identified as Hg-methylating. Highly abundant genera of Hg-methylators across all sampling sites included Geobacter and Syntrophales. Meanwhile, the phylum of Hg-methylators greatest in abundance across sampling sites was Desulfobacterota (Figure 6).

There was no significant difference in the average percentage of 16S rRNA ASVs that matched to Hg-methylating genera between Winooski River and Englesby Brook (p=0.3701). However, among all samples from both streams, the Riparian Land sites had significantly greater average percentage of 16S rRNA ASVs that matched to previously reported Hg-methylating genera (2.6%) compared to Riparian Water (1.0%) (p=0.03731), but not compared to Bay (2.4%) (Figure 7) [42]. While the negative control had some contaminants, there were only 13 ASVs found, as compared to an average of 1580 ASVs for the 15 samples from various sites (Figure 5), and no ASV contaminants found were identified as Hg-methylating genera.

3.3. Metagenome Analysis

Of the total 4,738,586 contigs assembled from our 6 metagenome samples, 5,426,668 genes were annotated. Of these, 46 genes (~0.001%) were identified as HgcA proteins using HMMER (v3.4) analysis (Table A3; Table A4; Table A5) [46]. These protein sequences were classified based on where they had been placed on the maximum likelihood reference HgcA tree. This analysis revealed that half of the sequences could only be identified to a superkingdom level (Table A5). Most of the assigned taxonomy of these proteins were Bacteria, PVC group, and “unclassified Deltaproteobacteria (miscellaneous).” Genera identified from four of the sequences were Clostridium, Syntrophus, Geobacter, and Methanoregula. The classifications of these sequences did not fully match the classifications identified through BLAST.

BLAST search of the identified protein sequences revealed that only 15 out of these 46 protein sequences had a top hit matching an existing HgcA protein sequence within the NCBI database (Table A4). Winooski RW and Winooski Bay had the greatest number of contigs matching to HgcA taxonomy. However, there was only one contig with 100% identity to the HgcA taxonomy identified. Some of the top hits for the protein sequences were mercury methylation corrinoid protein HgcA. However, there were also several top hits described as acetyl-CoA decarbonylase/synthase complex subunit gamma, acetyl-CoA synthase, acetyl- CoA synthase subunit gamma, carbon monoxide dehydrogenase, and “hypothetical protein.” The most abundant taxon identified as having HgcA proteins based on the metagenome analysis were Vicinamibacteria bacterium, with four contigs across the sampling sites being identified as matching to that species. Following that, contigs identified as having hgcA genes also frequently matched to Coriobacteriia bacterium and Deltaproteobacteria bacterium.

4. Discussion

4.1. Mercury Concentrations in Winooski River and Englesby Brook

THg concentrations within sediment samples were generally low and there was no significant difference in THg concentrations between the Riparian Land, Riparian Water, and Bay sites (Figure 2). THg concentrations reported in soil samples from other studies, such as sites near historical point sources of pollution, were around 1000 ng/g [9,49]. THg concentrations found in Vermont sediment samples were also typically higher than those from our study, though such studies are not directly comparable as they did not take sediment samples [26,50]. A possible reason for the generally low THg concentrations is that we sampled during the autumn season, when temperature and stream discharge are the lowest in the year, which are factors predicted to influence THg concentrations [26,51]. It is also possible that the intense floods of July 2023 may have resulted in deposition of sediments with low THg concentrations if any of the flood sediments were derived from a source area with low THg [52].

Comparing the two streams, Englesby Brook sites had higher THg concentrations with a mean of 41.94 ng/g. These results are consistent with previous findings which reported Englesby Brook as having the highest median THg concentration among 18 Lake Champlain tributary sites, as well as the highest annual areal flux at an average of 9.3 µg/m² from 2001-2009 [53]. Englesby Brook is also situated in a more urban location, and so has increased exposure to industrial sources of Hg, as well as concrete and asphalt surfaces increasing runoff of Hg and other environmental contaminants into the brook. In 2011, Burlington, the city in which 84% of the watershed of Englesby Brook falls into, was categorized as one of the municipalities in Vermont regulated by the Environmental Protection Agency and under the National Pollutant Discharge Elimination System (NDPES) [27]. Finally, due to being a much smaller stream compared to Winooski River, Englesby Brook likely experiences water level fluctuations at a much greater extent than Winooski Brook. This might increase Hg methylation because during conditions of low water level, sulfide can be oxidized to sulfate, which is then mobilized by potentially mercury-methylating SRB during conditions of high-water level [21,22].

4.2. Carbon and Nitrogen Concentrations in Winooski River and Englesby Brook

Overall, average carbon and nitrogen concentrations in dried sediment samples were slightly lower than previously described [5,54]. In Winooski River, the mean carbon concentration in sediment samples (9.97 mg/g dry wt) was an order of magnitude lower than that in wetland sediments [5], while Englesby Brook’s mean carbon concentration (24.3 mg/g dry wt) was similar (Table A1; Table A2). Nitrogen concentrations across all samples averaged 0.85 mg/g, with many samples below the detection limit, aligning with previous river floodplain studies, such as the Thur River [54], where the maximum nitrogen content was 1.8 g/kg. Most nitrogen-deficient samples were from Winooski River, suggesting it is more nitrogen-limited than Englesby Brook, likely due to Winooski River’s sandier sediment composition, as supported by previous research which found that sandy sediments contained less nitrogen [54]. In contrast, Englesby Brook sediment had more silt and clay.

There was no significant difference in average carbon concentration between the streams, but significant differences were observed between Winooski Bay compared to both Winooski Riparian Land and Winooski Riparian Water (Table A1). This was likely due to the high carbon content in Winooski Bay, where debris like wood pieces were present in the sediment samples. This higher carbon content in Winooski Bay suggests potential for elevated mercury (Hg) levels, as indicated by the positive correlation between total Hg (THg) and carbon (Figure 2) and the greater number of Hg-methylating contigs identified in the bay (Table A3 and A4). This significant positive relationship is consistent with previous studies showing a correlation between organic carbon and MeHg in sediments [49,55,56]. This relationship may also be explained through bioaccumulation [1]. Prior research [9] observed a positive correlation between MeHg and total carbon in a sediment site with decomposed fiber accumulation, where microbial processes such as carbon-fixation and sulfide oxidation were likely involved in MeHg production. Although there were significant differences in carbon concentration between Winooski sampling sites, no significant differences in mean THg concentrations were observed. This suggests that future studies should increase sample size and consider sampling at different seasons, to further determine the relationship between carbon and THg concentrations [26,51].

While there was a positive relationship between THg and carbon concentration in Englesby Brook sediment samples, it was not significant (Figure 3B), and the slope was less steep than in Winooski River. This could suggest that most mercury in Englesby Brook sediment remains inorganic, with less conversion to methylmercury (MeHg) and bioaccumulation compared to Winooski River. This is consistent with prior work[26] that found Englesby Brook had one of the lowest methylmercury fractions in the Lake Champlain basin. Our findings also show that Hg-methylating genes were more abundant and the percentage of 16S rRNA ASVs matching Hg-methylating microbes were greater in Winooski River sediments (Table A3; Figure 6).

4.3. Microbial Mercury Methylating Potential in Winooski River and Englesby Brook

The microbial community in the sediment samples showed high diversity at both the phylum and class level, similar to previous studies on Hg-methylating microbes in soil [57,58], with Gammaproteobacteria and Desulfomonilia being the most abundant classes. Several Hg-methylating genera were identified in the sediment samples through 16S rRNA ASVs (Figure 6), consistent with findings in contaminated sediments, which also showed a diverse range of Hg-methylating microbes at the class or genus level. The relative abundance of Hg-methylating genera in our samples was moderately high, exceeding 0.1% in relative abundance, which is considered significant in diverse communities like soil [59]. The most abundant Hg-methylating phylum was Desulfobacterota (Figure 6), which aligns with similar studies on Hg-methylating microbes within soil sediments [5,32,60]. Meanwhile the most abundant Hg-methylating genera were Geobacter, previously also found to be capable of Hg-methylation, and Syntrophales (Figure 6) [61].

The NMDS plot showed no clear microbial community differences between sites within each river (Figure 4), but 16S rRNA analysis revealed that Riparian Land and Bay areas had the highest relative abundance of Hg-methylating microbes, with RL having significantly greater relative abundance than RW (p=0.03731) (Figure 7). Previous studies show a significant correlation between Hg-methylating communities in sediment and MeHg concentration [58], suggesting that the presence of these microbes leads to increased MeHg production. While RL samples in Englesby Brook had higher THg concentrations than RW and Bay, no such difference was observed in Winooski River, likely due to lower carbon concentrations in the sediment, which may hinder the detection of bioaccumulated MeHg (Figure 2; Table A1).

In contrast to the 16S rRNA results, metagenomic analysis showed that RW samples had the most amino acid sequences matching Hg-methylating proteins, not RL samples (Table A3). This discrepancy may be due to the pooling of sediment samples for sequencing, which reduced the granularity of microbial diversity from different sites. The large number of sequences only classified at the superkingdom level when placing sequences on the maximum likelihood reference HgcA tree suggests that the HgcA proteins identified from HMMER search are highly divergent from known HgcA sequences. Meanwhile, the sequences matching to an unclassified Deltaproteobacteria aligns with our current knowledge about Hg-methylating microbes being dominated by sulfate-reducing bacteria [16,57,58].

The differing taxonomic classifications generated through placement of sequences onto a reference tree as compared to a BLAST search might lie in the differing algorithms used to obtain taxonomy. While pplacer aligns sequences of our samples to an existing model only limited to Hg-methylating taxa, BLAST uses pairwise alignment to compare two sequences [62]. Even top hits for BLAST searches were not very high, ranging from 69-100%, with only one sequence having 100% sequence identity (Table A4; Table A5). Further, phylogenetic placement onto a reference tree can be more accurate compared to similarity-based methods when dealing with sequences from metagenome samples with no exact or close matches within the existing database [63]. Phylogenetic placement can also account for greater variability of genes due evolutionary history and events such as horizontal gene transfer, gene duplication, etc. [64].

One surprising finding from the BLAST search of identified HgcA sequences was that the greatest number of sequences matched Vicinamibacteria, a genus recently identified as Hg-methylating [65]. This was not one of the genera detected during the 16S rRNA analysis using originally reported HgcA orthologs [14]. Vicinamibacteria, typically aerobic and neutrophilic [66], have been found in diverse environments, including peat fens [67], cadmium-polluted fields [68], and wastewater treatment plants [59]. Additionally, the abundant Deltaproteobacteria in our samples was expected, as this genus is well-known for Hg-methylation [5,58,69]. The wide distribution of hgcAB sequences supports the hypothesis of horizontal gene transfer of these genes among microbes [32]. This could be tested by analyzing the surrounding ATCG composition to determine if the genes are located on plasmids or areas susceptible to conjugation, transformation, or transduction.

Although only hgcA genes were analyzed to identify Hg-methylating microbes in sediment samples, as the presence of hgcA strongly correlates with Hg-methylating activity [58], hgcB genes are also essential for Hg-methylation. A positive correlation between hgcA abundance and MeHg content has been observed [70], and MeHg levels are significantly linked to sulfate-reducing bacteria (SRB) in other environments [69]. However, the presence of hgcAB genes doesn’t always indicate active Hg-methylation, as some microbes, like Acidobacteria, may have hgcA but not express it [32], highlighting the need for confirmation of methylating activity in a future study. Previous studies have also shown that the hgcAB genes are not crucial for survival or growth of the microbes capable of methylating, and mercury methylation cannot occur without the presence of mercury [14]. Additionally, several metagenomic contigs matched proteins linked to the acetyl-CoA pathway, likely due to the overlap between the Hg-methylation and acetyl-CoA pathways, as the corrinoid protein involved in Hg-methylation is also used in acetyl-CoA metabolism [14,17,32].

5. Conclusions

This study aimed to understand the presence of Hg-methylating microbes and key nutrients in sediment samples impacted by legacy mercury and sulfate deposition. By characterizing microbial communities, searching for Hg-methylating genes in sediment samples from Winooski River and Englesby Brook, Vermont, and analyzing the concentrations of total mercury (THg), carbon, and nitrogen, we gained insights into the relationship between mercury methylation and sediment nutrients. Notably, we discovered a positive correlation between carbon concentration and total mercury concentration in the Winooski River. Winooski River and Englesby Brook had a differentiated microbial community, with Winooski River having a higher relative abundance of Hg-methylating microbes based on 16S rRNA analysis. Additionally, while 16S rRNA analysis suggested that Riparian Land and Bay areas had a greater relative abundance of Hg-methylating microbes, metagenomic analysis indicated that Riparian Water samples had a greater number of sequences identified as hgcA genes. Metagenomics analysis also showed divergent Hg-methylating genes from taxa that were not previously categorized as Hg-methylating from the 16S rRNA analysis. Our study reveals the need for further investigation into mercury-methylating microbes as well as their functional capabilities to better understand these populations and their role in mercury cycling. These findings contribute to the growing understanding of mercury pollution in freshwater ecosystems and provide guidance for future studies on microbial mercury methylation.