Benchmarking IRcyc-A against other tools
We benchmarked IRcyc-A utilized in the Galaxy Europe server environment (Fig. 2E) against the existing metagenomics platforms MG-RAST23 and MGX24 that both, similar to IRcyc-A, provide a code-free graphic user interface without the necessity for user-provided access to supercomputing hubs. For the comparative analysis, we used a set of six publicly available soil metagenomes from US agroecosystems. The functional annotation against the SEED database in MG-RAST is capable of assigning reads to Fe(III)-reduction genes based on functional homology to the decaheme cytochrome and porin cluster of genes (MtrA-F) of Shewanella oneidensis. MGX on the other hand has the capacity to functionally annotate metagenomic reads against the TIGRFAM domain database25 thus also assigning reads to Fe(III)-reduction genes based on homology to the conserved domains TIGR03507 (OmcA/MtrC family decaheme c-type cytochrome), TIGR03508 (DmsE family decaheme c-type cytochrome including MtrA and MtrD), and TIGR03509 (MtrB/PioB family decaheme-associated outer membrane protein). Since both the MG-RAST-SEED database and TIGRFAM database used in MGX lack entries to Fe(II)-oxidation genes, we only compared the potential of the three databases to functionally annotate short reads to Fe(III)-reduction genes. Reads annotated by MGX as TIGR03508 (DmsE/MtrA decaheme) domain resulted in ~ 42% false-positive hits that did not match any of the so-far discovered MtrA genes provided by the exhaustive list found in the IRcyc-A database. The great number of false-positive hits found in functional annotation of metagenomic reads against the TIGRFAM domain database in MGX most likely results from a significant number of reads that belong to paralogous genes not included in the original TIGRFAM list of domains. Annotation against the SEED database in MG-RAST resulted in a small number of false positive hits – only 2% in line with their smaller number of matches. False-positive hits from the MGX and MG-RAST outputs were excluded for further downstream analyses.
Our findings show that IRcyc-A discovered 13-times greater number of reads mapping to the decaheme protein MtrA/DmsE/MtoA (1650 reads) compared to MG-RAST-SEED (132 reads), and 1.6-times more relative to the MGX-TIGRFAM (1059 reads) algorithm (Fig. 2B). That pattern was also apparent for other iron cycling genes (Fig. 2B). The greater discovery potential of IRcyc-A results from the more exhaustive list of full-length iron cycling protein gene sequences available within the IRcyc-A database. For instance, the IRcyc-A database contains MtrA/DmsE/MtoA entries from 1435 bacterial genomes compared to just one full length protein sequence, that of Shewanella oneidensis, in the SEED database of MG-RAST and just 13 used in the original seed alignment behind the conserved structure of the TIGR03508 domain.
Next, we performed taxonomic analyses of the MtrA reads reported by the three algorithms and found that in general the greatest number of taxa contributing to the MtrA pool in these soil metagenomes was found in searches using IRcyc-A followed by TIGRFAM-MGX, with the least found in SEED-MG-RAST searches (Fig. 2C). Phylum- and class-level taxonomic profiling revealed similar patterns between outputs of IRcyc-A and MGX-TIGRFAM with those of MG-RAST-SEED being most divergent, with heavy bias towards Proteobacteria, particularly Gammaproteobacteria (Fig. 2D). This bias is consistent with the fact that the SEED database version used in MG-RAST only contains MtrA-F sequences of the gamma-proteobacterium Shewanella oneidensis. This stands to demonstrate that the size and diversity of the database of iron cycling genes used in alignment searches ultimately determines the reliability and accuracy of downstream analyses.
Functional validation: local patterns in iron cycling predicted by IRcyc-A
• Acid sulfide soils in Finland
Acid sulfide soils are characterised by high amount of pyrite (FeS2) minerals found either as primary mineral inclusions in the underlying lithology or formed as a result of diagenesis of the soil substrate during coastal flooding and its subsequent reduced conditions26. They cover an estimated area of 17 million hectares, mainly in coastal regions of Australia, North America, Europe, and Asia27,28. The pyrite present in lower horizons (Fig. 3A) is unreactive under anoxic conditions. However, in contact with the oxygen atmosphere permeating upper soil horizons, pyrite interacts with oxygen in an oxidizing weathering process to form sulfuric acid and release ferrous iron (Eq. 1).
The resulting highly acidic conditions allow for iron to remain in the ferrous Fe(II) form and not become oxidized to Fe(III) despite the more oxic environment of the upper soil horizons29,30. Under these conditions, biological iron Fe(II)-oxidation by acidophilic bacteria6 and not abiotic oxidation is the dominant process converting Fe(II) to Fe(III) species30.
The acid sulfate soils at the Risöfladan experimental site in Finland are the result of a diagenesis process likely resulting from the rise in sea levels during the last 7500 years of the current Holocene era and subsequent rise of iron sulfide rich sediments accumulated in the anoxic sea zone to raise up to 100 m above sea level31. Inspection of the soil profiles at the site32 indicates the presence of three distinct types of sub-soil. The uppermost layer of the sub-soil contains an already oxidized zone (lowest pH, low ferrous iron) followed by an actively oxidizing zone (acidic pH slightly greater than that of the oxidized zone, highest measured ferrous iron levels) and lastly finishing in an unoxidized zone containing pyrite parent material (highest pH, lowest ferrous iron, Fig. 3A). To check if these recorded patterns in soil chemistry were aligned to changes in microbial iron cycling, we analysed the publicly deposited but as of yet unanalysed metagenomic and metatranscriptomic datasets associated with this experimental site32. Using the IRcyc-A database tool, we find that reads annotated to Fe(II) oxidation genes in the metagenomes predominantly mapped to the cytochrome c gene cyc2 utilized by iron oxidizers33,34. The total summed relative abundance of Fe(II) oxidation genes was greatest in the metagenomes of the oxidized and the actively oxidizing zones, whereas their abundance significantly dropped in the alkaline, unoxidized zone (Fig. 3A). Metatranscriptomic data available from this site contains information on the overall pattern of gene expression across the entire microbial community. Its analysis through the IRcyc-A database tool shows that the expression of Fe(II)-oxidation genes peaks in the actively oxidizing zone in line with the greatest levels of ferrous iron necessary for the lithoautotrophic metabolism of iron-oxidizing bacteria. Unlike their abundance in metagenomes, the relative expression of iron oxidation genes was relatively low in the metatranscriptomes of the uppermost oxidized zone. This apparent contrast between the two types of omics data likely results from the history of the oxidized zone. Iron oxidation was high in this zone at some point in the past. The subsequent exhaustion of ferrous iron substrate for oxidation consequently led to the cessation of the process of microbial iron oxidation but nevertheless left a transcriptionally inactive but metagenomically present potential for iron oxidation. The expression of Fe(III)-reduction genes was greatest in the actively oxidizing zone. This apparent synchrony between Fe(II)-oxidation and Fe(III)-reduction may be explained by the increasing concentrations of Fe(III) required by ferric reducers that is derived from the activity of iron oxidizers. This freshly generated Fe(III) will likely form more amorphous hydroxides and thus be more readily available for reductive dissolution by Fe(III)-reducers than the Fe(III) that had been oxidized in the past and that would already have achieved a certain level of crystallinity in its hydroxide minerals. Locally, the microsites enriched in Fe(III) with actively on-going pyrite and ferrous iron oxidation (see Eq. 1) would also consume great amounts of available oxygen perhaps creating microniches of anoxic conditions necessary for microbial Fe(III)-reduction. Interestingly, the expression of siderophore biosynthesis genes, in this case exemplified by the most abundant group in these microbiomes – that of the catecholate class of siderophores – was greatest in the metatranscriptomes of the bottommost unoxidized layer. These findings support the view that siderophore production and their iron-chelating activities are most favoured under conditions limiting the supply of iron available for immediate biological uptake by microorganisms – such as those of lowest levels of Fe(II) and alkaline pH as observed in the unoxidized zone.
• Intact versus drained peatlands in French Guiana
Peatlands cover an estimated 422 million hectares of land area worldwide35 and act as important sinks for carbon36. The peat found in peatlands contributes between a third and a half of the global reservoir of soil organic carbon – with temperate peatlands contributing 10–30%36 and tropical peatlands another ~ 20%37,38. Peat is a complex mixture of organic matter in different stages of decomposition formed under the anoxic conditions of these lands. Anoxic conditions develop owing the combined effects of high microbial respiration consuming oxygen and high soil water content limiting oxygen diffusion. Under these conditions, the processes of anaerobic microbial respiration including Fe(III)-reduction are favoured as alternative electron acceptors different to atmospheric oxygen (O2) have to be utilized39. The dissimilatory reduction of Fe(III)-containing minerals such as oxides goethite and ferrihydrite, phosphates such as vivianite and strengite, and carbonates such as siderite by microorganisms consumes electrons generated from the oxidation of electron donors – either organic (typically fermentation products such as lactate, acetate)40 or inorganic (hydrogen gas)39.
Draining is a method employed for converting peatlands for the purposes of agriculture, peat extraction, forestry, and human settlement41. Here we utilized a set of publicly available metagenomes from a study site in French Guiana41 to investigate the effects of peatland draining on the microbial Fe cycle. In terms of soil chemistry, the effects of draining resulted in predictable decreases in soil water content and soil organic carbon (Fig. 3B). The latter was likely caused by the increase of aerobic respiration, which allowed for faster decomposition of materials typically recalcitrant to decomposition under anoxic conditions. Using the IRcyc-A database, we reveal that draining was associated with a decrease in the relative abundance of Fe(II)-oxidation and Fe(III)-reduction genes in the microbial community (Fig. 3B) consistent with a more limited role for biogeochemical Fe cycling in response to the transition from anaerobic to aerobic metabolism. In contrast, transition from intact to drained peatlands was associated with a significant increase in the relative abundance of biosynthesis genes for all three major siderophore classes including the hydroxamate, catecholate, and mixed ligand class (Fig. 3B). These findings likely reflect the increased need of Fe(III)-chelating siderophores under more oxic conditions that would transition iron speciation from the highly mobile and directly available for microbial transport Fe(II) towards the Fe(III) forms intractable to direct biological utilization. Using IRcyc-A, we demonstrate that factors decreasing peatland soil water content such as draining but also probably climate change-driven drought would affect the biogeochemical cycling of iron and decrease the microbiome potential for iron oxidation and reduction while increasing the abundance of iron-chelating siderophore genes. These changes in the microbiome would likely be associated with a decline in the export of dissolved Fe to rivers from drained or dried peatlands.
• Biological iron oxidation as a mechanism driving silicate weathering and pedogenesis in tropical forest ecosystems in Australia
The Atherton Tablelands in Northeast Australia is a plateau developed on deep basaltic soils, the natural climax vegetation of which is luxuriant tropical rainforests. Basalt bedrocks in this system are geologically young, deposited mainly during the Quaternary period42. Mineralogically, these basalts are iron-rich42, with iron mostly in a ferrous form [Fe(II)] and contributed mainly by olivine and clinopyroxenes43. The oxidation of basaltic ferrous iron to ferric iron during weathering could contribute to pedogenesis and the substantial amount of secondary iron oxyhydroxide (mostly hematite44) minerals and acidity in these tropical soils. This is exemplified by the dissolution equations of the silicate iron-rich end-member of the olivine group – fayalite:
To investigate the role of biological iron oxidation of basalt as an acidogenic factor contributing to rock weathering and soil formation (pedogenesis) in this system, we deposited mesh bags containing fresh unweathered basalt grains (basalt bags, n = 16) or unweathered olivine (olivine bags, n = 17) in the top 0–10 cm of the secondary forest soils and left them to weather45 (Fig. 3C). After a year of in situ weathering, we recovered the bags and collected soil samples from the immediate surroundings of the deposited bags (n = 2). Next, we extracted total DNA from the collected weathered rock grains and soils and performed shotgun metagenome sequencing45. Here, we run these publicly available metagenomic sequencing data through the IRcyc-A pipeline. Our analyses show that both weathered basalt and soils exhibited significantly greater relative abundance of the Fe(II)-oxidation cyc2 genes than the microbial communities of weathered olivine (Fig. 3C).
Interestingly, both weathered basalt and olivine samples exhibited clear evidence for alteration, colour change, and development of iron oxide-rich rind (surface layer – see Fig. 3C). Comparing the pH of the three substrates in our study, it is clear that the weathering of basalt increased pH from 4.2 (surrounding soil) to 4.8 (inside basalt bags) but not to the same extent as olivine at pH 6.6 (inside olivine bags). These findings show that, in agreement with geochemical theory, the pH, but not total ferrous iron content (assuming majority as FeO and not Fe2O3 in fresh rock, basalt = 6.8% vs. olivine = 10.3% wt.) determined the extent of biological iron oxidation which is favoured by low pH29.
In soils, however, both acidity and total Fe are appropriate for biological iron oxidation (Fig. 3C), but most of the iron is present in ferric form as hematite. Consequently, the high abundance of cyc2 genes in these soils possibly reflect the redox fluctuations experienced during the wet season46. After prolonged wet weather, the water content of these tropical soils becomes conducive to microbial iron reduction of Fe(III) to Fe(II)46. During oxic conditions the thus accumulated Fe(II) can become biologically re-oxidized to Fe(III) by cyc2-carrying bacteria. However, given the large amounts of ferrous iron needed to support autotrophic growth by iron oxidation47, it seems likely that most iron oxidizers in highly weathered tropical forest soils (where most iron is ferric) are mixotrophs or heterotrophs that only partially complement their energy needs by cyc2-enabled iron oxidation – in an analogy to carboxydotrophic bacteria in soil48. In further support of that, the acidobacterial genus Ca. Acidoferrum, recovered as a MAG from soil in Panamanian tropical forests11, does contain cyc2 and cyc1 homologues but no carbon fixation pathway.
In contrast, during the weathering of basalt, there is no necessity for ferrous iron to be supplied by microbial iron reduction since most of the iron released from weathering is ferrous. At the pH of 4.8 (inside basalt rock bags), ferrous iron is relatively stable29 and locally abundant therefore allowing for biological iron oxidation by autotrophs49. Our generated data allowed us to examine this hypothesis by isolating several metagenome-assembled genomes (MAGs) from weathered basalt bags recovered from these tropical forest soils. One of these MAGs was identified to belong to a new species within the GJ-E10 lineage – a lineage of unnamed Burkholderiales bacteria previously shown to carry out iron oxidation under acidic conditions50. Genomic analysis of our GJ-E10-related MAG revealed the presence of homologues to cyc2 and cyc1 of Acidithiobacillus ferrooxydans and a full set of genes required for the Calvin cycle of autotrophic carbon dioxide fixation (Fig. 3C). Here, we propose a new taxonomic nomenclature for this generically named iron-oxidizing autotrophic genus GJ-E10 – Ca. Ferrobasaltibacterium, gen. nov. and for our newly discovered species – Ca. Ferrobasaltibacterium australiensis, sp. nov.
Taxonomic analysis of the universal single copy marker gene rpoA from these microbiomes through the IRcyc-A database show that Ca. Ferrobasaltibacterium was 10-fold more abundant in weathered basalt communities relative to those of soils and weathered olivine (Fig. 3C). These results strongly indicate that biological iron oxidation of basalt in these tropical soils requires basalt specialist autotrophic microbial assemblages and that iron oxidizing bacteria are important factor for silicate weathering and soil development under tropical conditions.
Global biome patterns in iron cycle genes
We selected soil microbiomes from a wide range of publicly available datasets with the final collection containing 193 soil shotgun metagenomes and 29 soil metatranscriptomes, representative of all major terrestrial biomes (Fig. 4A). A sub-collection of soil metagenomes (n = 15) representing specific treatment conditions were separated and the remaining metagenomes exemplifying untreated soil conditions were used for across-biome comparisons of iron cycling gene abundance.
The abundance of genes involved in the process of Fe(III)-reduction was highest (ANOVA, P < 0.001, post-hoc multiple comparisons) in soil metagenomes of peatlands and flooded grasslands (Fig. 4B). On average, the abundance of Fe(III)-reduction genes in the soil metagenomes of peatlands and flooded grasslands (mean = 0.79, SE = 0.04) was ~ 3.0-fold greater than the unweighted average across biomes (mean = 0.27, SE = 0.05). To establish the robustness of this observation, we compared the peatland metagenomes of several locations across the globe against the global average. Metagenomes of peatlands located in the United Kingdom, USA, Sweden, and French Guiana all exhibited means greater than the global unweighted average across biomes. Peatlands may thus represent a hotspot for biological Fe(III)-reduction consistent with their role in exporting larger amounts of soluble iron compared to non-peaty catchments (Fig. 1). From a metabolic perspective, organic matter oxidation coupled to the reduction of Fe(III) is most energetically favourable when alternative electron acceptors (O2, NO3−) are exhausted upon prolonged growth under anaerobic, low-redox conditions created by the high soil water content limiting oxygen diffusion and the high demand for O2 under the organic-rich peat. In line with these, the abundance of Fe(III)-reduction genes in soil metagenomes were lowest in desert and arctic deserts where available soil water and organic matter are at their lowest.
The abundance of genes involved in Fe(II)-oxidation was greatest (ANOVA, P < 0.001, post-hoc multiple comparisons) in soil metagenomes of acid iron sulfide soils (Fig. 4C). On average, the metagenomes of acid iron sulfide soils (mean = 0.74, SE = 0.15) exhibited 4.6-fold greater abundance in Fe(II)-oxidizing genes relative to the global across-biome average (unweighted mean = 0.16, SE = 0.06). This observation held for comparisons made between the two sites of acid iron sulfide soils – one in China and one in Finland – against the global average. Both of these acid sulfide soil sites exhibited mean Fe(II)-oxidation gene abundance 2.5-6.6-fold greater than the mean abundance in metagenomes across biomes. Iron(II)-oxidation is most commonly associated with acidic conditions under which Fe(II) remains stable even under highly oxic conditions. Iron sulfide soils thus provide the ideal conditions for biological Fe(II)-oxidation as the dissolution of their contained pyrite produces ferrous iron and the acidity (Eq. 1) necessary to maintain Fe(II) stable in oxic environments. In the absence of pyrite, acidophilic Fe(II)-oxidizers will be dependent upon the external acidity of the substrate (e.g., basalt weathered in tropical soils as in Fig. 3C) and natural fluctuations in soil redox conditions that (e.g., highly weathered oxisols in the tropics – Fig. 3C) generate Fe(II) from soil Fe(III)-oxyhydroxides through microbial Fe(III)-reducing activity. These conditions are best exemplified by the soils of tropical rainforests, peatlands, and tundra where high soil acidity and/or fluctuations in soil redox are commonplace. Indeed, our global biome analyses confirm that soil metagenomes from these three biomes exhibited some of the greatest abundance of Fe(II)-oxidizing genes after that hereby reported for acid sulfide soils (Fig. 4C).
Next, we analysed the across-biome patterns in the metagenomic gene abundance of different classes of siderophore biosynthesis genes. Genes involved in the biosynthesis of catecholate siderophores (e.g., enterobactin, vibriobactin, bacillibactin etc.) were greatest in the soil microbiomes of tropical grasslands (savannas), acid iron sulfide soils, peatlands, and tropical forests, while agroecosystems and temperate grasslands showed lower abundances (Fig. 4D). Lower yet were the catecholate gene abundances of soil microbiomes of tundra, deserts, and arctic deserts. These results may highlight the dual role of catecholates serving as Fe(III)-chelating agents produced in response to iron deficiency51,52 as well as anti-oxidant molecules protecting against reactive oxygen species (ROS)53. ROS are commonly produced in soil as a result of the respiratory activity of soil microorganisms54. Biotic ROS undergo a catalytic reaction with soil Fe(II) and Fe(III)-(oxyhydr)oxide minerals that generate the highly reactive OH• radicals.
These Fenton and Fenton-like reactions are particularly enhanced in soils experiencing high levels of microbial activity, fluctuations in redox state, Fe(II) production via microbial Fe(III)-reduction, and are rich in poorly-crystalline amorphous iron minerals54. These findings help explain the high abundance of genes involved in the production of ROS-scavenging catecholate siderophores in microbiomes typically associated with high biological availability of iron – e.g., peatlands, acid iron sulfide soils, tropical forest soils etc. or large fluctuations in redox state such savannas and to a lesser extent tropical forests experiencing alternating dry and wet seasons (Fig. 4D). These provide a novel insight between Fe(II)-oxidizers and Fe(III)-reducers on the one hand affecting the levels of reactive iron species and Fenton-based ROS generation54 and catecholate siderophore producers on the other hand that can inactivate the produced ROS species53. In addition, the Fenton-based chemistry of iron cycling and the role of catecholates in these cases will ultimately determine the amount of ROS species that can impact soil organic matter stability and decomposition55. These findings help us gain further insight into the coupling between the global Fe and C cycling.
The metagenomic relative abundance of genes implicated in the production of hydroxamate siderophores (e.g., desferrioxamine B and E, ferrichrome, aerobactin, malleobactin56 etc.) was highest in agroecosystems and deserts closely followed by temperate forest, temperate grasslands, tropical grasslands, and tropical forests (Fig. 4E). The lowest abundance in hydroxamate biosynthesis genes was observed in the soil metagenomes of acid iron sulfide soils, peatlands, tundra, and arctic deserts. Genes involved in the synthesis of mixed ligand siderophores (e.g., pyoverdine, mycobactin, yersiniabactin etc.) followed a similar pattern to that of hydroxamate siderophores in that agroecosystem soil metagenomes had among the greatest abundance, whereas lowest abundances were recorded for the soil microbiomes of arctic deserts, deserts and tundra. Furthermore, we performed correlational analyses that revealed a cluster of strong positive correlations between the relative abundance of Fe(III)-reducing, Fe(II)-oxidation, and catecholate biosynthesis genes across biomes (Fig. 4F). A second cluster of positive correlations was formed between the relative abundance of genes implicated in the biosynthesis of hydroxamate, mixed ligand, and carboxylate siderophores (Fig. 4F). The two clusters were negatively correlating with each other. Our findings can be understood as a result of microbiomes being actively shaped by iron availability and ROS species. For instance, catecholate siderophores, Fe(II)-oxidation genes, and Fe(III) reduction genes are highest in abundance in soil microbiomes subjected to high biological availability of iron in soil allowing for extensive biogeochemical cycling of iron and ROS species evolution. This is consistent with previously suggested views that the processes of iron reduction and oxidation frequently occur simultaneously or cyclically7. In contrast, the genes for the majority of siderophore classes (hydroxamates, mixed ligands, carboxylates) peak in metagenomic abundance where Fe availability is more limited as indicated by their negative association with Fe(II)-oxidation and Fe(III)-reduction genes. These findings are in line with a major role of siderophores in ameliorating Fe deficiency57.
Taxonomic profiling of target Fe cycling gene groups across biomes (Fig. 5) suggest that members of the Acidobacteriota, many of which uncultured58, dominate the gene pools involved in biogeochemical cycling of iron (oxidation and reduction; Fig. 5A,B) in topsoils, whereas Actinobacteriota dominate the production of hydroxamate siderophores (Fig. 5C). Actinobacteria and Acidobacteriota together with the group of Proteobacteria are among the most important phyla in terms of iron cycling and chelation (Fig. 5A-C). Actinobacteria and Acidobacteria are also among the most abundant bacterial phyla in soil microbiomes across the globe, making up to ~ 60% of all bacteria in soil (Fig. 5D, panel 1). We observed a strong negative correlation between the abundance of these groups in soil microbiomes across biomes (Fig. 5D, panel 1). Furthermore, as expected based on their contribution to the pools of Fe cycling genes, the abundance of Acidobacteriota in the metagenome correlated positively with the abundance of Fe(III)-reduction and Fe(II)-oxidation genes (Fig. 5D, panels 2 and 3). Similarly, the abundance of Actinobacteriota correlated positively with the abundance of hydroxamate siderophore biosynthesis genes in the metagenome (Fig. 5D, panel 4). These findings suggest that the abundance of these two keystone phyla in soil microbiomes is associated with iron availability, chelation and redox cycling, pointing to iron as an overlooked factor in shaping soil microbial communities.