Apr 21, 2024
Drought re
Nature Microbiology volume 8, pages 1480–1494 (2023)Cite this article 1878 Accesses 28 Altmetric Metrics details Drought impacts on microbial activity can alter soil carbon fate and lead to the loss
Nature Microbiology volume 8, pages 1480–1494 (2023)Cite this article
1878 Accesses
28 Altmetric
Metrics details
Drought impacts on microbial activity can alter soil carbon fate and lead to the loss of stored carbon to the atmosphere as CO2 and volatile organic compounds (VOCs). Here we examined drought impacts on carbon allocation by soil microbes in the Biosphere 2 artificial tropical rainforest by tracking 13C from position-specific 13C-pyruvate into CO2 and VOCs in parallel with multi-omics. During drought, efflux of 13C-enriched acetate, acetone and C4H6O2 (diacetyl) increased. These changes represent increased production and buildup of intermediate metabolites driven by decreased carbon cycling efficiency. Simultaneously,13C-CO2 efflux decreased, driven by a decrease in microbial activity. However, the microbial carbon allocation to energy gain relative to biosynthesis was unchanged, signifying maintained energy demand for biosynthesis of VOCs and other drought-stress-induced pathways. Overall, while carbon loss to the atmosphere via CO2 decreased during drought, carbon loss via efflux of VOCs increased, indicating microbially induced shifts in soil carbon fate.
Microorganisms regulate terrestrial carbon (C) cycling in fundamental ways1, including transforming soil C into gaseous compounds that can escape to the atmosphere, primarily as CO2 via microbial heterotrophic respiration. However, microbes also produce volatile organic compounds (VOCs) as metabolic intermediates, signalling molecules and secondary metabolites2,3. In fact, volatile metabolites represent an often overlooked subset of the complete soil metabolome4,5, and although their emissions to the atmosphere only represent a small soil C loss, they contribute substantially to atmospheric chemistry including ozone formation and cloud condensation nuclei6. Therefore, characterizing microbe-mediated C flow along the soil–atmosphere continuum is critical for understanding the fate of soil C and VOCs under projected environmental changes, including drought.
Drought stress induces well-characterized microbial physiological responses that impact C metabolism, such as biosynthesis of protective molecules (for example, osmolytes and extracellular polymeric substances) to preserve cell integrity7,8 and concentrate resources9,10. Production of these biomolecules is energy intensive and may divert resources from biomass synthesis9, leading to decreased growth-fuelling heterotrophic respiration and CO2 emissions11,12. Drought also induces changes in soil C composition and availability13,14,15, further impacting microbial activity16,17, for example, by inducing shifts in substrate utilization18. Overall, it remains unclear how drought-induced shifts in microbial metabolism and soil C composition influence C allocation to volatile metabolites, which can mitigate drought stress in plants19. Moreover, soil water content has a strong impact on VOC emissions from soils20, including tropical soils21,22, perhaps driven by drought impacts on microbial VOC biosynthesis and/or consumption. Characterizing changes in microbial C cycling and allocation is particularly important in tropical rainforest soils where droughts will probably occur more frequently and last longer due to climate change23,24.
Detecting shifts in microbial C cycling and allocation within complex metabolic networks encompassing competing production and consumption pathways is challenging. This complexity can be navigated by tracking microbe-mediated C flow through soils using isotopically labelled central metabolites. Position-specific 13C-glucose and/or 13C-pyruvate labelling has been used to track microbial C allocation to CO2 and biomass in soil mesocosms25,26,27. However, studies on drought in the field and allocation to VOCs are missing28. The direct metabolic information derived from isotope labelling can be contextualized using powerful, information-rich constraints provided by ‘omics approaches that profile the gene content, gene expression and metabolomes of soil microbiomes26. Together, these approaches can uncover metabolic drivers of shifting microbial C cycling and allocation in soils under drought.
In this study, we performed a comprehensive assessment of how drought impacts soil microbial C allocation to CO2 and VOCs in an artificial enclosed tropical rainforest at Biosphere 2 using an integrated approach, combining position-specific 13C-pyruvate labelling with metatranscriptomics, metagenomics and metabolomics. Biosphere 2 (Oracle, Arizona) is a 12,700 m2 steel and glass-enclosed building harbouring five distinct biomes, including a 1,900 m2 tropical rainforest—an enclosed ecosystem with variable topography and microhabitats where rain and temperature can be controlled. In 2019, a 65-day drought was imposed on the tropical rainforest to study ecosystem-scale impacts of dry-down conditions29. With the assumption that the C1 position of pyruvate would enter CO2 via decarboxylation during biosynthesis or respiration (tricarboxylic acid (TCA) cycle), and the C2 position would enter VOCs (as biosynthesis) or CO2 during respiration26,27, we aimed to infer microbial C allocation to biosynthesis vs TCA cycle during both ambient and drought conditions. We hypothesized that during drought: (1) microbial C allocation to the TCA cycle would decrease as C allocation is diverted to increased biosynthesis of stress compounds, including VOCs and (2) changes in C allocation would be associated with shifts in metabolic composition and gene expression.
We tracked 13C from the first and second positions of pyruvate (13C1-pyruvate and 13C2-pyruvate, respectively) into CO2 and VOCs (Fig. 1) within soil chambers located among three sites across the tropical rainforest at Biosphere 2 (Extended Data Fig. 1). During drought, total soil CO2 efflux and its 13C-enrichment (δ13CCO2) decreased after injection of 13C1-pyruvate and 13C2-pyruvate relative to pre-drought, with higher δ13CCO2 from 13C1-pyruvate compared with 13C2-pyruvate chambers during both pre-drought and drought (Fig. 2a and Extended Data Fig. 2a). As expected, this pattern was reflected in the cumulative 13C-CO2 effluxes (Fig. 2b), which decreased significantly during drought by 54.3% from chambers receiving 13C1-pyruvate (t-value = −5.85, d.f. = 7, P = 0, linear mixed effect (LME)) and 47.6% from chambers receiving 13C2-pyruvate (t-value = −3.92, d.f. = 7, P = 0.0058, LME).
Framework depicting how first (C1; blue) and second (C2; pink) carbon (C) positions of pyruvate from 13C1-pyruvate and 13C2-pyruvate, respectively, are expected to be allocated into CO2 and biosynthetic (including VOCs) pathways. As indicated by the arrow thickness, we expect most C1 to be released as CO2 during biosynthesis and TCA cycle, and C2 to be either released as CO2 in the TCA cycle or used for biosynthesis of biomass and products. The difference in 13C-CO2-C1 − 13C-CO2-C2 approximates the total C allocation to biosynthesis, while the ratio 13C-CO2-C1 / (13C-CO2-C1 + 13C-CO2-C2) approximates the proportion of internal C allocation to biosynthesis. *Some 13C from 13C1-pyruvate may end up in the TCA cycle due to anaplerotic CO2 assimilation, which would lead to a slight overestimate in the proportion of C allocation to biosynthesis and overall C allocation to biosynthesis. Figure adapted with permission from refs. 26,27, Elsevier.
a, δ13C of CO2 efflux over time from 12 h pre to 48 h post pyruvate injection from both 13C1- and 13C2-pyruvate. b, Cumulative 13C-CO2 soil efflux from 0 to 48 h post pyruvate injection. c, Overall C allocation to biosynthesis (calculated as 13C-CO2-C1 − 13C-CO2-C2). Continuous emission data were binned to 3 or 6 h intervals. d–f, Fluxes of 13C-enriched (13C/(12C + 13C)) acetic acid (d), acetone (e) and C4H6O2 (diacetyl+) (f) from 0 to 48 h post pyruvate injection. g–i, Cumulative flux of 13C-acetic acid (g), 13C-acetone (h) and 13C-diacetyl+ (i) from chambers that received 13C1- and 13C2-pyruvate. j, Pie charts of the percentage of 13C from pyruvate allocated to biosynthesis or the TCA cycle (or other which could represent unmetabolized pyruvate, CO2 fixation or uncharacterized pathways) for pre-drought (P) and drought (D). For all time series (a,d–f), each point represents a single measurement per soil chamber (~64 measurements), of which there were 3 replicates each that were injected with either 13C1- or 13C2-pyruvate for each site (total n = 18 each for pre-drought and drought), and lines show the smoothed data with the surrounding shaded area showing ±s.e.m. Boxplots represent Q1–Q3, centre lines indicate the median and whiskers extend to the minimum and maximum values, exclusive of outliers (black points), of cumulative flux (measured per soil chamber) across all sites and timepoints (b,g–i) (n = 54 each for pre-drought and drought) or overall C allocation to biosynthesis (measured per set of C1/C2 chambers) across sites and specific timepoints as indicated (c) (n = 9 sets each for pre-drought and drought). P values in b, c, g–i are based on linear mixed dffects models. *P < 0.05, **P < 0.01, ***P < 0.001.
Source data
Contrary to our hypothesis, microbes maintained their allocation to energy production during drought-induced stress. While overall C allocation to biosynthesis was higher during pre-drought conditions between 3 and 18 h post pyruvate addition (Fig. 2c), internal partitioning of C to biosynthesis vs TCA cycle did not change (Extended Data Fig. 2b).
To further reveal C allocation strategies, we tracked 13C from pyruvate into 13C-enriched (13C/(12C + 13C)) volatile compounds and were only able to detect three, all of which are central metabolites: acetone, acetic acid and C4H6O2 (Fig. 2d–f). Cumulative effluxes of 13C-acetic acid and 13C-acetone increased significantly by factors of 2.6 and 25.3, respectively, during drought from chambers receiving 13C2-pyruvate (t-values = 2.83 and 4.30, respectively; d.f. = 7, P < 0.05, LME) (Fig. 2g,h), indicating synthesis along metabolic pathways after C1 decarboxylation to CO2. Emissions of acetone also increased during low moisture in a previous study in the tropical rainforest21, and these results suggest a role of microbes in the production of acetone and consequent emissions during drought. Acetic acid and C4H6O2 also showed 13C-enriched continuous efflux (Fig. 2e,f) from chambers that received 13C1-pyruvate, with cumulative effluxes that increased significantly by factors of 1.9 and 3.5, respectively, during drought (t-values = 3.45 and 2.76, respectively; d.f. = 7, P < 0.05, LME) (Fig. 2g). 13C-enrichment of acetic acid from chambers receiving 13C1-pyruvate indicates active acetogenesis, where acetogens fix two CO2, in this case 13C-CO2, and four H2 molecules to form acetyl-CoA, which is further reduced to acetate30 (Fig. 1), an intermediate of central metabolism that can be protonated in soils and volatilized as acetic acid (hereafter referred to as acetate). C4H6O2 probably comprised diacetyl and/or one other unidentifiable compound (henceforth referred to as diacetyl+) based on structural identification of compounds with this formula from nearby locations (Extended Data Fig. 3). 13C may have been allocated to other non-volatile metabolites; however, no such 13C-enrichment was detected by 1H-nuclear magnetic resonance (NMR).
Drought induced a readjustment of overall C allocation. The proportion of 13C from pyruvate allocated to biosynthesis dropped from 41.0 to 17.3% and to the TCA cycle from 21.1 to 11.1% during drought (Fig. 2j). Despite this overall decrease in C allocation to biosynthesis, 13C-enriched emissions of acetate, acetone and diacetyl+ represented a greater proportion of biosynthetic products during drought, increasing from 0.0040% to 0.15% (Fig. 2j). Meanwhile, 13C from pyruvate that did not track to emitted 13C-CO2 (‘other’), possibly representing unmetabolized pyruvate and/or autotrophic and heterotrophic CO2 fixation, increased from 37.9 to 71.6% during drought (Fig. 2j). Total CO2 fixation would be expected to decrease CO2 emissions by up to 5.6%, depending on soil type30,31, suggesting the potential for only a small impact on our estimates of microbial C allocation to biosynthesis vs TCA cycle. Overall, these data suggest an increase in biosynthesis of VOCs during drought, despite an overall decrease in microbial C allocation to biosynthesis.
To characterize active microbial metabolic pathways that could cycle acetate, acetone and diacetyl+ in soils, we utilized a gene-centric approach using metatranscriptomics and metagenomics data to identify: (1) specific genes that may be driving the 13C-enriched VOC emissions (see Fig. 3a for possible metabolic pathways leading from pyruvate to acetate, acetone and diacetyl, and Extended Data Fig. 6 for associated gene expression (acetate and acetone only)) and (2) overall gene expression patterns that reflect ecosystem-scale microbial responses to drought. Drought and 13C-pyruvate injection (that is, comparison between 0 (before injection), 6 and 48 h post-injection timepoints) significantly affected overall gene expression profiles (P < 0.05, permutational multivariate analysis of variance (PERMANOVA)) (Extended Data Fig. 4a) and taxonomic composition of active microbial communities (Extended Data Fig. 4c–e) but had no impact on microbial function potential (Extended Data Fig. 4b) or taxonomic composition (Extended Data Fig. 4f–h). Only location (site) impacted microbial functional potential (P = 0.003, PERMANOVA) (Extended Data Fig. 4b). This demonstrates that fluctuating gene expression was driven by changes in microbial activity and not by drastic shifts in community composition. We point out the caveat that these results are based on known genes and pathways within the Kyoto Encylopedia of Genes and Genomes (KEGG) database, potentially leading to bias against some organisms that are underrepresented in the databases, such as protists32.
a, Pyruvate metabolism pathways that lead to 13C-enriched VOCs and CO2, shown with blue and pink circles, to represent incorporation of 13C-pyruvate from the C1 and/or C2 position, respectively. Colour code for compound names: grey names without brackets indicate compounds that were not measured, grey names with brackets indicate VOCs that were not emitted, plain black names indicate VOC fluxes that were emitted but not 13C-enriched and bold black names indicate VOC emissions that were 13C-enriched (acetate, acetone, diacetyl+). Red arrows represent a potential PDH bypass route of acetyl-CoA production from pyruvate to acetate to acetyl-CoA. Green and orange gene names signify down- or upregulation during drought, respectively, compared to pre-drought conditions. b, Changes in gene expression across all sites (for pre-drought: 0 and 48 h (n = 6), 6 h (n = 5); for drought: n = 6 for each timepoint) as log2 fold-change (log2FC) of drought vs pre-drought for acetate, acetone and diacetyl+ (C4H6O2) at 0, 6 and 48 h post pyruvate addition. Genes that encode enzymes involved in production are in black bars and consumption in grey bars. P values in b are based on DESeq2 and are FDR corrected. ^, anaplerotic CO2 assimilation via pyc (pyruvate carboxylase) leading to some 13C from 13C1-pyruvate entering the TCA cycle and being released as 13C-CO2. acxA–C, acetone carboxylase subunits A–C; ilvB, acetolactate synthase; cs, citrate synthase; ald1, NAD+-dependent secondary alcohol dehydrogenase; OAA, oxaloacetate; PEP, phosphoenolpyruvate. *P < 0.05, **P < 0.01, ***P < 0.001.
Source data
Both 13C-enriched acetate and acetone effluxes from chambers that received 13C2-pyruvate (acetate-C2 and acetone-C2, respectively) were driven by a combination of genes encoding for production and consumption, as well as soil moisture, which together explained 91% of acetate-C2 and 76% of acetone-C2 efflux (partial least square regression (PLSR)) (Table 1). Soil moisture was a major driver of both acetate-C2 and acetone-C2 effluxes (VIP = 1.02 and 1.00, w1 = −0.42 and −0.44; PLSR) (Table 1), indicating that emission of these microbially produced VOCs was partly driven by drought-induced moisture changes. Previous research found that decreased soil moisture is associated with increased movement of volatile compounds through the soil33.
Microbial gene expression that drove acetate-C2 efflux primarily included acetate-producing enzymes, three of which were upregulated during drought (Fig. 3b). Expression of ACH1 (K01067), encoding for an acetate-producing acetyl-CoA hydrolase which was upregulated during drought (P = 8.4 × 10−14, differential analysis using DESeq2) (Fig. 3b and Extended Data Fig. 5b), was the largest microbial driver of acetate-C2 efflux (VIP = 1.51; w1 = 0.51; PLSR) (Table 1). ACH1, only present in eukaryotes33,34,35, including fungi33,35, produces acetate and CoA from acetyl-CoA36, and can play a role in acetyl-CoA regulation. The two other largest microbial drivers of acetate-C2 efflux were expression of poxB (K00156; encoding pyruvate dehydrogenase (PDH)-quinone) and poxL (K00158; encoding pyruvate oxidase) (VIP = 0.90 and 0.92, w1 = 0.17 and 0.18, respectively; PLSR) (Table 1). The poxB gene, first identified in Escherichia coli37,38, encodes for the protein PDH-quinone, which catalyses the direct conversion of pyruvate to acetate and CO2 via oxidative decarboxylation, with quinone or menaquinone as the electron acceptor37,38. Quinone then shuttles the electrons to the electron transport chain where O2 is the final electron acceptor, thereby producing ATP39. The strong upregulation of poxB during drought suggests that this gene may improve microbial fitness during drought, perhaps to gain extra energy for biomolecule production, as discussed below.
Microbial gene expression that drove acetone-C2 efflux was limited to three genes: two encoding for acetone-consuming enzymes and one for acetone production. Expression of adc (K01574), encoding for the acetone-producing enzyme acetoacetate decarboxylase, was the largest microbial driver of acetone-C2 efflux (VIP = 1.05, w1 = 0.61; PLSR) (Table 1); however, adc expression did not change during drought (Fig. 3b). The gene adc is part of the acetone-butanol-ethanol (ABE) fermentation pathways in Clostridium acetybuticulum and related species40,41; however, we did not see the expected concurrent 13C-enriched ethanol or butanol emissions, indicating either the presence of other routes of acetone production, or maintained consumption of ethanol and butanol under drought conditions. The other two genes that drove acetone-C2 efflux were adh1 (K18382) and acmA (K18371), which encode for the acetone-consuming enzymes alcohol dehydrogenase (producing propanol, which was not detectable due to its extreme fragmentation during ionization with proton-transfer-reaction mass spectrometry (PTR–MS))42 and acetone monooxygenase (producing methyl acetate), respectively (VIP = 0.84 and 0.64, respectively; PLSR) (Table 1). AcmA was also downregulated during drought, along with three other genes encoding for acetone-consuming enzymes (acxA (K10855), acxB (K10854) and acxC (K10586)) (P < 0.05; DESeq2) (Fig. 3b). While C3H6O2 was detected, possibly indicating methyl acetate, it was not 13C-enriched. It is possible that methyl acetate was further broken down into methanol and 13C-acetate via methyl acetate hydrolase (acm; K18372); however, acm was not a major driver of acetone-C2 or acetate-C2 efflux. Therefore, the overall significant decrease in acetone-C2 efflux was driven by both increased production and decreased consumption during drought, signifying that changes in VOC fluxes from soils depend not only on production but also on internal consumption.
For diacetyl+-C2 efflux, gene expression and soil moisture accounted for only 23% of total variation (Table 1). Diacetyl is formed when ilvB (K01652), encoding for the enzyme acetolactate synthase, catalyses the conversion of pyruvate to CO2 and acetolactate, which is further converted non-enzymatically to diacetyl43. During drought, ilvB was downregulated (P = 0.029; DESeq2) (Fig. 3b) and its expression was not a major driver of diacetyl+-C2 efflux. Only the expression of butA (K03366) was a microbial driver of diacetyl+-C2 efflux (VIP = 1.05; PLSR) (Table 1). ButA encodes for meso-butanediol dehydrogenase, an enzyme that plays a diacetyl-consuming role producing acetoin (Fig. 3a), although there was no significant change in expression of butA during drought (Fig. 3b). The unidentifiable C4H6O2 compound that was present along with diacetyl (Extended Data Fig. 3) may have also driven 13C-enriched emissions, therefore making it difficult to assess genes involved in diacetyl+ cycling.
To place 13C-enriched emissions of acetate, acetone and diacetyl+ and the metabolic pathways involved in the cycling of these volatile metabolites into the context of overall microbial responses to drought and substrate availability, we further characterized bulk soil metabolites and microbial gene expression.
Metabolite composition shifted between pre-drought and drought conditions, which separated along the first component (PC1), explaining 44.7 and 62.3% of variation in small (<200–300 Da) compounds measured using 1H-NMR (Extended Data Fig. 7a) and relatively larger compounds (>200 Da) characterized using Fourier-transform-ion-cyclotron-resonance MS (FTICR–MS) (Extended Data Fig. 7b), respectively. Small compounds, mainly representing primary metabolites, that decreased during drought included several amino and keto acids (oxiosocaproate, t-value = −5.09; alanine, t-value = −3.72; formate, t-value = −6.25; glycine, t-value = −3.12; leucine, t-value = −3.16; phenylalanine, t-value = −2.71; pyroglutamate, t-value = −3.09; pyruvate, t-value = −5.31; uracil, t-value = −2.94; and valine, t-value = −3.34; d.f. = 7, P < 0.05, LME), while only trehalose, a common osmolyte produced by bacteria during times of water stress7, increased during drought (t-value = 3.00; d.f. = 7, P = 0.020, LME) (Supplementary Table 1). Drought can induce increased soil concentrations of trehalose14,15, yet the impact on amino acids is variable, with some studies finding an increase13,14,15 and some a decrease44. While we did not detect non-volatile 13C-enriched primary metabolites with 1H-NMR, decreased concentrations of amino acids during drought could indicate decreased C allocation to biosynthesis of biomass or enhanced degradation of proteins.
For relatively larger compounds, representing microbial metabolism and substrate availability, recalcitrant lignin-, condensed hydrocarbon- and tannin-like compounds increased (t-values = 4.17, 4.16, 2.43, respectively; d.f. = 35, P < 0.05, LME) and bioavailable protein- and amino-sugar-like compounds, along with lipid-like compounds, decreased during drought (t-values = −5.37, −3.92, −5.44, respectively; d.f. = 35, P < 0.05, LME) (Extended Data Fig. 7b,c). One possible explanation for this pattern is that microbes preferentially consumed bioavailable compounds during drought, while there was a concurrent decrease in labile C replenishment from root exudates from drought-stressed plants29.
Next, we looked at overall shifts in microbial metabolic pathways of co-expressed gene modules using weighted gene co-expression network analysis (WGCNA). We identified a total of nine modules, four of which were associated with condition (pre-drought or drought) and acetate-C2 and acetone-C2 effluxes (Extended Data Fig. 8a,b).
The pre-drought-associated red module was enriched in C-cycling pathways. The association of the red module with pre-drought conditions was demonstrated by its negative correlation with (r = −0.6, P = 5.0 × 10−4, Pearson correlation coefficient (PCC)) (Extended Data Fig. 8b) and downregulation during (t-value = −4.24, d.f. = 30, P = 3.8 × 10−4, LME) (Fig. 4a) drought. The enriched central C metabolic pathways interwoven within the red module included C metabolism, butanoate metabolism (includes acetone and diacetyl cycling), and pyruvate and methane metabolism (both include acetate cycling) (Fig. 4a). These overlapping central C metabolic pathways suggest efficient use of intermediate C compounds produced, including volatile intermediates such as acetone and acetate, and are supported by negative correlations with acetate-C2 and acetone-C2 efflux (r = −0.32 (P = 0.090, not significant (NS)) and −0.47 (P = 0.010), respectively; PCC) (Extended Data Fig. 8b). Therefore, efficient C cycling in the red module led to rapid consumption of metabolic intermediates, preventing losses to the atmosphere under moist conditions.
a–d, Top: eigengene expression at 0, 6 and 48 h post pyruvate addition during pre-drought (0 and 48 h (n = 6), 6 h (n = 5)) and drought (6 and 48 h (n = 6), 0 h (n = 5)) conditions for the subset of modules: carbon-efficient, red (a); stress-response, pink (b); fungal, green (c) and archaeal, magenta (d). Expression values are arbitrary units. Each point represents one sample; boxes represent Q1–Q3 with centre line indicating median, and bars extend to maximum and minimum values, excluding outliers. Bottom: metabolic network of enriched KEGG pathways within the indicated module. Central nodes represent the pathway and each branch represents the KO group within that pathway. P values in a–d are from linear mixed effects models between pre-drought and drought (as indicated after ‘Condition’) or 6 or 48 h vs 0 h (as indicated by lines between timepoints); *P < 0.05, **P < 0.01, ***P < 0.001.
Source data
The drought-associated pink module played a role in osmotic-stress adaptation, quinone production and acetate/acetone biosynthesis. The pink module was positively correlated with (r = 0.8, P = 1.0 × 10−7, PCC) (Extended Data Fig. 8b) and upregulated during (t-value = 7.59, d.f. = 30, P < 1.8 × 10−8, LME) (Fig. 4b) drought, and it was also potentially involved in acetate and acetone production, as indicated by its positive correlations with acetate-C2 and acetone-C2 efflux (r = 0.29 (P = 0.010, NS) and r = 0.43 (P = 0.020), PCC) (Extended Data Fig. 8b) and inclusion of the acetate-producing poxB and poxL genes. Furthermore, the pink module was enriched in starch and sucrose metabolism; alanine, aspartate and glutamine metabolism; and ribosome pathways (Fig. 4b). Biosynthesis of trehalose within the starch and sucrose metabolism pathway corresponds to the increase in trehalose concentrations observed during drought (Extended Data Table 1). At the sub-pathway level, biosynthesis of ubiquinone, a type of quinone, was also enriched, indicating a possible association with osmotic stress. We hypothesize two explanations that connect quinone biosynthesis, poxB expression and acetate biosynthesis as stress-response mechanisms during drought. First, a link between the electron transport chain and osmotic regulation is possible. As an immediate microbial response to osmotic stress, microbes actively pump in K+, which then promotes trehalose biosynthesis45. A supercomplex of H+-ATPase and K+ pumps may form in the cellular membrane as part of the electron transport chain46. Therefore, increased production of quinone and transport of electrons to the electron transport chain could be linked with K+ influx, leading to further osmotic adjustment. Second, due to the high demand for protective biomolecule synthesis (that is, trehalose), there is a high demand for energy production during drought, which may be limited by lower TCA cycle activity (Fig. 2j). To address this need, poxB may act as a bacterial PDH bypass route for production of acetyl-CoA (in combination with acs)47, the precursor to many downstream secondary metabolites48 (Fig. 3a). While experiments in Corynebacterium found that the PDH bypass route was not essential for growth38, poxB contributes to the aerobic growth efficiency in glucose-limited conditions in E. coli47. It is, therefore, possible that under times of stress or within different bacterial species, this route improves fitness by increasing acetyl-CoA production to meet energy demands.
The drought-associated green and magenta modules were enriched in fungal and archaeal metabolic pathways, respectively. Both modules were significantly correlated with (r = 0.54 and 0.36, respectively; P ≤ 0.05, PCC) (Extended Data Fig. 8b) and upregulated during (t-value = 3.94 and 2.26, respectively; d.f. = 30, P < 0.05, LME) (Fig. 4c,d) drought; however, while the green module was positively correlated with acetate-C2 and acetone-C2 (r = 0.62 and 0.77, respectively; P < 0.001, PCC) (Extended Data Fig. 8b) and contained the acetate-cycling genes ACH1 and adc indicating its association with acetate and acetone efflux, the magenta module showed no correlation with 13C-enriched VOC emissions or genes of interest. The abundance of fungal genes within the green module’s enriched KEGG pathways indicates a potential role for fungi in the drought response; specifically, the green module was enriched in genes encoding for eukaryotic ribosomal subunits within the ribosome pathway, as well as eukaryotic-specific enzymes (NADH:ubiquinone oxidoreductase and F-type ATPase) within the oxidative phosphorylation pathway (Fig. 4c) (KEGG pathway database49). Therefore, the association between the green module and acetate and acetone efflux, as well as fungal-specific genes, suggests a potential role of fungi in acetate and acetone fermentation during drought. While previous studies have found that fungal communities are generally more resistant to stresses, such as drought, compared with bacteria50,51,52, here we did not see a significant increase in total active fungi during drought based on taxonomic classification (0.16 ± 0.7 and 1.1 ± 0.3% in relative abundance during pre-drought and drought conditions, respectively; t-value = 0.63, d.f. = 136, P = 0.52, LME); however, specific fungal taxa shifted in activity during drought, with Ascomycota increasing on average from 44.3 to 68.6% of the fungal community during drought (t-value = 3.02, d.f. = 31, P = 0.0050, LME) (Extended Data Fig. 4e), with Fusarium comprising a majority of this phylum particularly during drought at 48 h post pyruvate addition (52%). Specific studies examining fungal production of acetone are very limited; however, acetone production by Penicillium brevicompactum has been detected53, and Fusarium spp. have been found to produce several volatiles, including acetone and acetate54. In contrast, the magenta module was enriched in pathways specific to archaea. For example, the enriched DNA replication and RNA polymerase pathways (Fig. 4d) included polymerases specific to archaea. Despite the association between the drought-associated magenta module and drought, there was no concurrent increase in active archaeal abundance during drought (t-value = 0.87, d.f. = 276, P = 0.37, LME) (Extended Data Fig. 4d). Collectively, these results suggest that microbial acetate and acetone biosynthesis during drought may be associated more with fungal rather than archaeal metabolic activity, and we suggest further studies to experimentally test this hypothesis.
Here we show how microbially mediated soil C allocation strategies shift from typical (pre-drought) to drought conditions using a combination of 13C-pyruvate labelling and multi-omics in a unique controlled-drought experiment performed within an enclosed tropical rainforest. We found a distinct subset of emitted volatiles affected by drought that were formed from pyruvate, hence our findings are representative of pyruvate-dependent volatile emissions from forest floor areas that are not covered by understory vegetation. Under typical pre-drought conditions, there was a balance of C allocated to biosynthesis vs energy, and we were able to track 13C into several volatile compounds, all representing primary metabolic intermediates (acetate, acetone and diacetyl+). We hypothesize that most of the C was allocated to the biosynthesis of biomass as well, because the 13C-enriched volatile intermediates only made up a tiny fraction of the total C allocated to biosynthesis. Therefore, production of these volatiles was matched with immediate recapturing of these molecules during efficient C cycling and utilization. During the stressed conditions imposed by drought, there was a disruption to C allocation pathways. While the balance between C allocated to biosynthesis vs energy remained constant compared to pre-drought conditions, biosynthesis shifted (from presumably biomass) to stress biomolecules, including trehalose and quinone. Furthermore, metabolic pathway shifts led to higher release of ‘leaky’ central metabolites (indicated by increased emissions of 13C-enriched acetate, acetone and diacetyl+), enhanced by C allocation to VOC production (that is, acetate) and reduced capacity to consume and recoup VOC-C (that is, acetone)55. Notably, this decreased C cycling efficiency during drought could be amplified by increased air-filled pore space, isolating microbial communities to microhabitats within the soil matrix that are cut off from an aqueous flow of non-volatile metabolites and substrates12. This signifies that during drought, microbes switch from investing in more-stable carbon pools such as biomass9 to stress compounds that are prone to be quickly utilized by microbes upon rain rewet56. Overall, microbial survival responses to drought shift soil C cycling by interrupting C storage pathways, increasing C allocation to VOCs and loss in efficiency of recouping C from volatile intermediates, thus indicating potential shifts in soil carbon fate.
Biosphere 2 in Oracle, Arizona is a 12,700 m2 steel and glass-enclosed building harbouring five distinct biomes, including a 1,900 m2 tropical rainforest—an enclosed ecosystem with variable topography and microhabitats with the ability to control rain and temperature inside, providing an optimal setting for studying drought effects57,58,59. This artificial tropical rainforest harbours approximately 70 species of trees and shrubs forming a canopy and understory, with soils that represent biogeochemical cycles present in natural rainforests60. In late 2019, a 65-day drought experiment was conducted as part of the Water, Atmosphere and Life Dynamics (WALD) campaign29. During ambient (pre-drought) conditions, rainfall events were simulated by spraying 15,000 l of irrigation water from the top of the Biosphere 2 tropical rainforest at a frequency of 3 d a week, with the last rainfall event before the drought occurring on 7 October 2019. During the drought, ambient temperatures were maintained between 20 and 26.7 °C in the lowland region. For more detailed characterization of the WALD drought experiment, please see ref. 29.
Bulk soil was labelled with position-specific (C1 or C2) 13C-pyruvate (henceforth referred to as 13C1-pyruvate and 13C2-pyruvate, respectively) to track C allocation into CO2 and VOCs—a method adapted from one used in plants61. Similar 13C-pyruvate labelling was performed in other regions of the tropical rainforest during the WALD campaign, including plant roots62 and leaves63. Three replicates each of 13C1-pyruvate or 13C2-pyruvate labelling were performed per site (Site 1, Site 2 and Site 3; n = 6 per site), representing a vast expanse of the lowland region (Extended Data Fig. 1a) within soil chambers (Extended Data Fig. 1b) before (12–19 September) and during (7–19 November) drought. The night before labelling, two automatic chambers were placed onto the corresponding soil collars, which contained no plants but might have had small amounts of fine roots and were covered with a rain-out shelter. Each morning, these collars were labelled at around 10:00 by placing a 5 × 5 cm stencil with 1 × 1 cm openings into one side of each chamber (placed on different sides of the chamber during pre-drought and drought) and adding 100 μl of 13C1-pyruvate or 13C2-pyruvate solution (40 mg ml−1) (Cambridge Isotopes, CLM-8077-PK and CLM-8849-PK) to each opening to a depth of 1 cm (25 injections), for a total of 100 mg (Extended Data Fig. 1b). The stencil was then removed before soil chamber gas measurements.
Soil moisture and temperature were measured using a portable probe and LabQuest viewer (Vernier). For the pre-drought condition, soil moisture and temperature data were collected on 1 and 9 October and for drought on 11 and 18 November for a subset of collars (P11, P21, P23, P32 and P33 (Site # (P1, P2 or P3), replicate # (1, 2 or 3)), except for 18 November when all collars from the experiment were tested (see Extended Data Fig. 1a). Soil moisture measured near the labelling sites decreased from 26.0 ± 6.9 to 13.8 ± 2.6% (P < 0.001) between pre-drought and drought conditions, respectively. Soil temperature did not change significantly and ranged from 23.0 ± 0.6 during pre-drought to 23.3 ± 1.2 °C during drought.
Before labelling, collars were measured at typical temporal resolution (~2 h) overnight. To capture any rapid changes in gas fluxes due to changes in microbial activity after pyruvate addition, measurement intervals were increased to high frequency (30 min) directly before labelling. After gas fluxes were expected to equilibrate, approximately 8 h post pyruvate labelling (~18:00), measurement intervals were decreased to low frequency (50 min) and remained at this frequency until measurements were stopped at 48 h post labelling, resulting in ~64 measurements for each soil collar. For each measurement, the automatic chamber closed over the collar for a total of 10 min (pre-purge, 2.5 min; measurements, 6.5 min; post-purge, 1 min). Fluxes were measured using an automated multiplexed Licor soil flux system (Licor 8100, Li-8150 16-port multiplexer and Lic 8100-104 long-term chambers with opaque lids, Licor). The system was coupled to a Picarro G2201-i analyser (Picarro) to measure CO2 and isotopic composition and a proton-transfer-reaction time-of-flight mass spectrometer (PTR–TOF 8000, Ionicon) for VOCs (including 13C-VOCs). The PTR–ToF–MS sampled the sub-flow from the soil flux system at 30 standard cubic centimeters per minute (sccm). Perfluoroalkoxy tubing was used for the soil flux system, for the subsampling line and for the PTR inlet, with the aim to minimize the release and the retention of the VOCs from and on the tubing surface64. The PTR settings were as follows: inlet temperature was 60 °C, drift voltage was 600 V, drift temperature was 60 °C and drift pressure was 2.2 mbar, resulting in an E/N (E being the electric field strength and N the sample gas number density) ratio of 137 Td (Townsend). The time resolution was 10 ms and the mass range was up to 500 amu. The PTR–ToF was operated in the H3O+ mode, therefore only compounds having proton affinity higher than water (697 kJ mol−1) underwent proton-transfer reactions and could be detected on their protonated mass to charge ratio (m/z), which includes most VOCs. Measured ions were attributed to chemical formulae and specific chemical species based on the exact protonated m/z. PTR–TOF data were processed using the software PTRwid65. To account for possible variations in reagent ion signals, measured ion intensities were normalized to the H3O+ counts in combination with the water-cluster ion counts66. At midnight, automatic calibrations were performed using standard gas cylinders containing different multi-VOC component calibration mixtures in Ultra-High Purity (UHP) nitrogen (Apel-Riemer Environmental). For a detailed description of the calibration setup, see ref. 29. The concentrations of compounds included in the standard were calculated with an uncertainty of ≤23%. Concentrations of compounds not included in the calibration standard cylinders were calculated by applying the kinetic theory of proton-transfer reaction67,68 with an uncertainty of ≤50%.
Select additional soil experiments were conducted with a Vocus PTR–ToF (TOFWERK)69 coupled to a custom-built gas chromatograph (GC) (Aerodyne Research)70. The GC contained an integrated two-stage adsorbent-based thermal desorption preconcentration system for in situ collection of VOCs before separation on the chromatographic column. For preconcentration, a multibed (Tenax TA/Graphitized Carbon/Carboxen 1000, Markes International) sorbent tube was used for the first stage of sample collection; this tube was then subjected to a post-collection water purge before the sample was thermally desorbed to a multibed focusing trap (Tenax, Carbopack X, Carboxen 1003, Markes International) before injection onto the GC column. For this study, the GC was equipped with a 30-m Rxi-624 column (Restek, 0.25 mm i.d., 1.4 µm film thickness) which resolves non- to mid-polarity VOCs in the C5–C12 volatility range before ionization in the PTR detector. The GC–PTR sampled from in situ soil gas probes71 on an alternating timed schedule. The GC–PTR can speciate structural isomers and help identify some unknowns by matching to calibrated retention times.
CO2 fluxes and their isotopic composition were calculated on the basis of data from the Picarro isotope analyser. Fluxes were calculated with linear and exponential models fitted to each individual chamber measurement. A deadband of 30 s was used for each measurement to allow for mixing in the chamber. The linear models were calculated on the basis of the first 120 s and for the exponential model, the full closure period of 6.5 min was used. Fluxes were quality controlled visually. Fluxes based on the exponential fit were used preferentially but were replaced by the linear-fit flux where necessary. The isotopic composition was calculated on the basis of the individual efflux rates of the 12C-CO2 and 13C-CO2 isotopologues. These are reported separately by the gas analyser and linear fits based on the first 120 s were used to calculate efflux rates. Due to the high enrichment of the labelled soil CO2 efflux, this method was found to be more reliable compared with conventional Keeling-plot approaches. The isotopic composition of the CO2 efflux was then calculated from the ratios of 12C-CO2 to 13C-CO2 efflux rates and normalized to the Vienna Pee Dee Belemnite (VPDB) scale (δ13CCO2 = ((13C/12C)CO2)/((13C/12C)VPDB) – 1). C isotope fluxes were quality controlled visually for each individual chamber and outliers removed manually.
13C-CO2 emitted from chambers that received 13C1-pyruvate is formed as C is decarboxylated via pyruvate dehydrogenase (PDH) to form acetyl-CoA or via an alternate decarboxylation reaction leading to biosynthesis of compounds (for example, biomass, secondary metabolites, VOCs), while 13C-CO2 emitted from chambers that received 13C2-pyruvate is primarily formed during decarboxylation in the TCA cycle during energy production (Fig. 1) (modified from refs. 26,27). Using this concept, we can approximate both total and relative (proportion of total) C allocated to biosynthesis by calculating 13C-CO2-C1 − 13C-CO2-C2 and 13C-CO2-C1/(13C-CO2-C1 + 13C-CO2-C2), respectively, where 13C-CO2-C1 and 13C2-CO2-C2 are efflux from chambers that received 13C1-pyruvate or 13C1-pyruvate, respectively. Within each site, there were three sets of ‘C1’ and ‘C2’ chambers located next to each other (Extended Data Fig. 1), and the calculations above were made for each set. To facilitate calculation of total and relative C allocation over time, continuous flux data were binned and averaged across 0–3 (3 h), 3–6 (6 h), 6–12 (12 h), 12–18 (18 h), 18–24 (24 h), 24–30 (30 h), 30–36 (36 h), 36–42 (42 h) and 42–48 (48 h). These calculations are approximations due to two reasons. First, heterotrophic and autotrophic CO2 fixation, including anaplerotic reactions31,72, leads to 13C-CO2 being immediately re-consumed after production, causing an underestimate of actual 13C-CO2 emissions by up to 5.6% (based on soil CO2 fixation rates calculated from a forest soil31). Second, due to anaplerotic CO2 assimilation, some 13C from 13C1-pyruvate may have entered the TCA cycle and been released as 13C-CO2, causing an overestimate of biosynthesis. However, assuming constant CO2 fixation routes in pre-drought and drought conditions, this would not have impacted the conclusions we made from our calculations.
The isotopic composition of VOC flux rates was calculated by applying the linear model to calculate the rate of change in the fractional abundance of 13C (13C-VOC/(13C-VOC + 12C-VOC)). For each 6.5 min chamber measurement, a deadband of 30 s was used to allow for mixing in the chamber and the linear model was applied to the successive 10 data points collected at 10 s intervals. To identify pyruvate metabolism pathways, we focused on those VOCs that showed 13C isotopic enrichment in their emissions. To help clarify metabolic pathways and direction of reactions, we also considered the fluxes of VOCs independent of their 13C-enrichment and metabolite concentrations immediately up- or downstream of the 13C-enriched VOCs (Fig. 3a). For a complete analysis of all VOCs that showed changes in flux patterns across the drought experiment, see ref. 73.
Cumulative 13C-CO2 and 13C-VOC effluxes from 0 to 48 h post 13C-pyruvate injection were calculated using the area under the curve ‘auc’ function within the FLUX package (v.0.3)74,75 in R on the basis of the total amount of 13C-pyruvate added, we could determine what percentage of total 13C ended up in CO2 or VOCs.
For metagenomics, metatranscriptomics and metabolomics, soil samples were collected just before 13C-pyruvate labelling (0 h) then at 6 and 48 h post labelling. For 0 h collection, soil (~ 8 g) was collected directly outside of the stencil using a 2.25 cm diameter sterile metal ring pushed into the soil to 2 cm depth and placed into a sterile 50 ml tube. For 6 and 48 h collections, soil samples were collected using the same method as for 0 h, but inside the metal frame where pyruvate labelling occurred. After each soil collection, samples were immediately brought to the lab and allocated for different downstream analyses: 1 g stored at −20 °C for metabolomics and 2 g in Lifeguard soil preservation solution (Qiagen, 12868–1000) stored at −80 °C for DNA/RNA extractions.
RNA and DNA were co-extracted from 1 g of soil using the RNeasy Powersoil Total RNA kit (Qiagen, 12866–25) coupled with the RNeasy Powersoil DNA Elution kit (Qiagen, 12867–25) following the manufacturer’s protocol and eluted in 100 μl of kit-provided elution buffer. RNA and DNA concentrations and quality were measured using a Qubit 4 fluorometer (Thermo Fisher) and NanoDrop spectrophotometer (Thermo Fisher). RNA was further treated with DNAse (DNAse Max, Qiagen, 15200–50) to remove any DNA contamination. Total RNA and DNA were sent to the Joint Genome Institute (JGI; Berkeley, California) for library prep and sequencing.
Water extractions were performed on samples for NMR, followed by solid phase extraction (SPE) for FTICR–MS. Water extraction procedures for bulk metabolite characterization were followed76,77,78. Briefly, 1 g of soil and 5 ml of double deionized water were vortexed in a 15 ml centrifuge tube, sonicated for 2 h at 21 °C, then centrifuged for 20 min. One ml of water-extractable organic C (WEOC) was removed and stored at −80 °C for NMR analysis at the Pacific Northwest National Laboratory (PNNL), while 4 ml of WEOC was saved for downstream SPE.
SPE is an extraction technique used to clean and concentrate samples for the isolation and analysis of organic compounds79. Four ml of WEOC were acidified to pH 2 with 1 M HCl, then passed through methanol (MeOH)-preactivated Bond Elut PPL (Priority PolLutant, Agilent, 12255002) barrel cartridges that contain a macroporous styrene-divinylbenzene crosslinked polymer to retain polar organic compounds, with a vacuum set no higher than −5 psi. Next, the cartridges were washed three times with 9 ml 0.01 M HCl to remove impurities, air dried with a pressure air hose for 2–3 min and finally rinsed with 1.5 ml MeOH in a slow dropwise flow rate into 2 ml vials. The eluate was stored at −80 °C and sent to PNNL for FTICR–MS analysis.
1H-NMR bulk metabolite characterization was performed on WEOC. Samples (180 µl) were combined with 2,2-dimethyl-2-silapentane-5-sulfonate-d6 (DSS-d6) in D2O (20 µl, 5 mM) and thoroughly mixed before transfer to 3 mm NMR tubes. Resonances corresponding to 13C labelling were identified by visual inspection, comparing labelled spectra to unlabelled spectra. Once differences were identified, the molecular location and quantification of 13C incorporation were determined by J-coupling pattern analysis and the ‘linefitting’ tool in MestReNova (v.14.2.014)80, respectively; however, 13C-labelled metabolites could not be identified in our samples. NMR spectra were acquired on a Bruker Advance III spectrometer operating at a field strength of 17.6 T (1H ν0 of 750.24 MHz) and equipped with a 5 mm Bruker TCI/CP HCN (inverse) cryoprobe with Z-gradient and at a regulated temperature of 298 K. The one-dimensional (1D) 1H spectra were acquired using a Nuclear Overhauser Effect Spectroscopy (NOESY) pulse sequence (noesypr1d). The 90° H pulse was calibrated before the measurement of each sample with a spectral width of 12 ppm and 1,024 transients. The NOESY mixing time was 100 ms and the acquisition time was 4 s, followed by a relaxation delay of 1.5 s during which presaturation of the water signal was applied. Time domain free induction decays (72,114 total points) were zero filled to 131,072 total points before Fourier transform, followed by exponential multiplication (0.3 Hz line-broadening) and semi-automatic multipoint smooth segments baseline correction. Chemical shifts were referenced to the 11H methyl or 13C signal in DSS-d6 at 0 ppm. The 1D 11H-NMR spectra of all samples were processed, assigned and analysed using Chenomx NMR suite 9.2 (Chenomx) with quantification of spectral intensities of compounds in the Chenomx library relative to the internal standard. Candidate metabolites present in each of the complex mixtures were determined by matching chemical shift, J-coupling and intensity information of the experimental signals against signals of the standard metabolites in the Chenomx library. Signal to noise ratios (S/N) were measured using MestReNova with the limits of quantification and detection equal to an S/N of 10 and 3, respectively. Standard 2D experiments such as 11H/13C - heteronuclear correlation (HSQC) or 2D 11H/11H Total Correlation spectroscopy (TOCSY) further aided corroboration of several metabolite identifications where there was sufficient S/N.
A 12-Tesla Bruker FTICR mass spectrometer was used to collect high-resolution mass spectra of SPE-filtered WEOC by direct injection for secondary metabolite characterization. Approximately 100 µl of water extract was mixed with MeOH (1:2) before injection to enhance ionization. Samples were directly injected into a standard Bruker electrospray ionization (ESI) source. The instrument settings were optimized by tuning on a Suwannee River fulvic acid standard (International Humic Substances Society), and the instrument was flushed between samples using a mixture of water and MeOH. Blanks (HPLC grade MeOH) were analysed at the beginning and end of the day to monitor potential carry over from one sample to another. The ion accumulation time varied to account for differences in C concentration between samples. For each sample, 144 individual scans were averaged and internally calibrated using an organic matter homologous series separated by 14 Da (CH2 groups). The mass measurement accuracy was <1 ppm for singly charged ions across a broad m/z range (100–1,000 m/z). The mass resolution was ∼240 K at 341 m/z. The transient was 0.8 s. Data Analysis software (BrukerDaltonik v.4.2) was used to convert raw spectra to a list of m/z values, applying the FTICR–MS peak picker module with an S/N threshold set to 7 and absolute intensity threshold set to the default value of 100. Putative chemical formulae were then assigned using Formularity software81 on the basis of the following criteria: S/N > 7, mass measurement error <1 ppm and taking into consideration the presence of C, H, O, N, S and P and excluding other elements82. To ensure consistent formula assignment and eliminate mass shifts that could impact formula assignment, all sample peak lists were aligned to each other. The following rules were implemented to further ensure consistent formula assignment: (1) picking formulae with the lowest error between predicted and observed m/z and the lowest number of heteroatoms and (2) requiring the presence of at least four oxygen atoms for the assignment of one phosphorus atom82. The chemical character of thousands of peaks in each sample’s ESI FTICR–MS spectrum was evaluated using van Krevelen diagrams, with biochemical compound classes reported as relative abundance values on the basis of counts of C, H and O as follows: lipids (0 < O:C ≤ 0.3 and 1.5 ≤ H:C ≤ 2.5), unsaturated hydrocarbons (0 ≤ O:C ≤ 0.125 and 0.8 ≤ H:C < 2.5), proteins (0.3 < O:C ≤ 0.55 and 1.5 ≤ H:C ≤ 2.3), amino sugars (0.55 < O:C ≤ 0.7 and 1.5 ≤ H:C ≤ 2.2), lignin (0.125 < O:C ≤ 0.65 and 0.8 ≤ H:C < 1.5), tannins (0.65 < O:C ≤ 1.1 and 0.8 ≤ H:C < 1.5) and condensed hydrocarbons (aromatics; 0 ≤ 200 O:C ≤ 0.95 and 0.2 ≤ H:C < 0.8)78. Analysis of FTICR data was performed using MetaboDirect (v.0.2.7)83 to create profiles of biochemical compound classes and principal component analysis (PCA) plots.
Metagenome and metatranscriptome libraries were created at the JGI. Plate-based DNA library preparation for Illumina sequencing was performed on the PerkinElmer Sciclone NGS robotic liquid handling system using the Kapa-HyperPrep library preparation kit (Kapa Biosystems). Sample genomic DNA (200 ng) was sheared to 300–500 bp using a Covaris LE220 focused-ultrasonicator. The sheared DNA fragments were size selected by double-SPRI, then the selected fragments were end-repaired, A-tailed and ligated with Illumina-compatible sequencing adaptors (from IDT) containing a unique molecular index barcode for each sample library.
For metatranscriptome libraries, ribosomal (r)RNA was removed from 100 ng of total RNA using Qiagen FastSelect 5S/16S/23S for bacterial rRNA depletion (and additional FastSelect plant and/or yeast rRNA depletion) (Qiagen) with RNA blocking oligo technology. The fragmented and rRNA-depleted RNA was reverse transcribed to create first-strand complementary (c)DNA using Illumina TruSeq Stranded mRNA Library prep kit (Illumina), followed by the second-strand cDNA synthesis which incorporated dUTP to quench the second strand during amplification. The double-stranded cDNA fragments were then A-tailed and ligated to JGI dual indexed Y-adapters, followed by an enrichment of the library by 10 cycles of PCR. Metatranscriptomics of sample P35SSC1_191108_c were not completed due to librarly prep and sequencing issues.
The prepared libraries were quantified using KAPA Biosystems’ next-generation sequencing library qPCR kit and run on a Roche LightCycler 480 real-time PCR instrument. Sequencing of the flowcell was performed on the Illumina NovaSeq sequencer using NovaSeq XP V1.5 reagent kits and S4 flowcell, following a 2 × 151 indexed run recipe.
Sequence filtering, assembly, mapping and annotation were performed at JGI84. BBDuk (v.38.94), included in the BBtools package85, was used to trim reads that contained adapter sequences, homopolymers of G’s of the size 5 or more at the ends of reads and reads where quality dropped to 0. BBDuk was also used to remove reads that contained 1 or more ‘N’ bases, had an average quality score across the read of <10, or had a minimum length of ≤51 bp or 33% of the full read length. Reads were mapped with BBMap (v.38.44), included in the BBtools package, to masked human, cat, dog and mouse references at 93% identity, and common microbial contaminants were removed for downstream analysis. Filtered reads were error corrected before assembly using bbcms (v.38.90), included in the BBtools package and assembled using metaSPAdes (v.3.15.2)86 with a minimum contig of 200 bp (for metatranscriptome libraries, there was an additional removal of rRNA during filtering and assembly was performed using MegaHit (v.1.2.9)87). Filtered reads were then mapped back to contigs using BBMap to obtain coverage information. See Extended Data Table 2 for sequencing depth and mapping statistics. Feature prediction was next performed on assembled contigs by using tRNAscan-SE (v.2.0.6)88 to predict transfer RNAs, INFERNAL (v.1.1.3)89 to identify non-coding RNAs and rRNAs, CRT-CLI90 to identify CRISPR regions, and Prodigal (v.2.6.3)91 and GeneMarkS-2 (v.1.07)92 to identify protein-coding genes (CDSs). Functional annotation was performed on CDSs with KEGG Orthology (KO) terms93. Taxonomy for each CDS was determined using the best LAST94 hits of CDSs, then contigs were classified on the basis of the majority rule (the lowest taxonomic rank that at least 50% of CDSs on the contig matched to84). Gene copies per KO were calculated as the average coverage of the contigs each gene was predicted from, multiplied by the number of genes in the KO84.
Active metabolic pathways were determined by mapping KOs to KEGG pathways using the KEGG pathway mapper tool93. VOC cycling genes were chosen on the basis of their being immediately up- or downstream from the VOC in the metabolic pathway, as well as being downstream from pyruvate. Gene symbols used in the paper are based on the KEGG database. Whether the gene was involved in production or consumption of VOC was assessed by examining the reaction on the KEGG database and on the basis of the literature.
All statistical analyses were performed in R (v.4.0.2)95. First, to compare 13C-CO2 and 13C-VOC (acetate, acetone and diacetyl+) cumulative fluxes across pre-drought and drought conditions, LME modelling was performed using the LME function in the NLME package (v.3.1)96 with ‘Site’ and ‘ID’ (location within each site) included as random variables to control for these confounding factors. Next, LME was used to compare differences in total and relative C allocation to biosynthesis (13C-CO2-C1–13C-CO2-C2 and 13C-CO2-C1/(13C-CO2-C1 + 13C-CO2-C2), respectively) between pre-drought and drought conditions. The following model equation was used for both of the above LME analyses (full scripts are included in GitHub at github.com/linneakh/SoilPyruvate):
WGCNA was performed on metatranscriptomic data to identify modules of co-expressed genes. First, gene copy data were normalized in DESeq2 (v.1.30.1)97 using the variance stabilization transformation (VST) method and samples were clustered to detect outliers, which removed only one sample (P15SSCI_191107_c). Next, the normalized data were analysed using the WGCNA package (v.1.7)98,99 in R employing the ‘blockwiseModules’ function with the following settings: soft-thresholding power of 6, minimum module size of 80, minimum KME of 0.35 and a ‘signed’ topology overlap matrix. Pearson correlation coefficients (r) were calculated between module eigengenes and the following: (1) condition (values for pre-drought and drought set to −1 and 1, respectively) and (2) 13C-enriched VOC fluxes (acetate-C2 and acetone-C2) averaged at 0–2, 5–7 and 46–48 h post pyruvate addition to correspond with metatranscriptome soil collection times of 0, 6 and 48 h, respectively.
To identify expressed genes that were contributing to 13C-enriched VOC efflux, PLSR analysis was performed using the pls package (v.2.8)100. This multivariate statistical test identifies predictive variables that contribute to a response variable while allowing for collinearity between variables. Three PLSR models were created using the 13C-enriched flux (acetate-C2, acetone-C2 and diacetyl+-C2) as a single response variable for each model. To match flux measurements to gene expression data, we averaged VOC flux measurements for 0–2, 5–7 and 46–48 h post pyruvate addition to correspond with metatranscriptome soil collection times of 0, 6 and 48 h, respectively. As predictive variables, the VST-normalized gene copy data were used from collars that received 13C2-pyruvate (Site2-P4, Site2-P6 and Site3-P4) from across three timepoints post pyruvate injection (0, 6 and 48 h) during pre-drought and drought conditions (n = 18). To account for the correlation between acetate-C2 and acetone-C2 at the 0, 6 and 48 h timepoints, which is probably due to the interconnection of the pathways for cycling of these compounds in metabolic pathways, acetate-C2 flux was added as a predictive variable to the acetone model, and acetone-C2 flux was added as a predictive variable to the acetate model. Furthermore, to improve the models, we included soil moisture (%) which likely impacts the movement of VOCs through the soil and emissions to the atmosphere8. The optimal number of components for each model was determined by ‘leave one out’ cross-validation, calculating the root mean squared error in cross-validation (RMSECV) using 0–10 components and finding the number of components that had the lowest RMSECV value. The degree of importance of each predictive variable was assessed using the variable importance in prediction (VIP) values, which were calculated from the loadings, weights and scores of each variable101 using the plsVarSel package (v.0.9.10)102. Additionally, to determine the significance of the model, cross-validated predictive residuals were compared to the residuals of the null model (that is, using the mean of the response variable) using an F-test. By comparing F-values for each model and their corresponding degrees of freedom, a P value could be calculated to determine whether the model was significantly different from the null. To determine which predictive variables significantly contributed to each model, the predictor variable with the lowest VIP was removed and a new model was created. This was repeated until the model with the highest % of variation explained (R2 × 100) was obtained, indicating the predictive variables that were the largest drivers of the response variables. The variables removed from each model were as follows: acmB (K18372) for the acetate model; acmB (K18372), acxC (K10856), acxA (K10855) and acxB (K10854) for the acetone model; and ilvB (K01652) and BDH (K0004) for the diacetyl+ model. The optimal numbers of components for each model were determined as follows: 6 components for acetate (R2 = 0.51, 0.19, 0.16, 0.03, 0.03 and 0.002 (total R2 = 0.91)), 1 component for acetone (R2 = 0.76) and 2 components for diacetyl (R2 = 0.15 and 0.08 (total R2 = 0.23)), with corresponding RMSECV values of 42.6, 62.5 and 31.6, respectively. This resulted in the following P values: 0.50 for acetate, 0.005 for acetone and 0.93 for diacetyl+ models.
To find differences in metabolomic composition and metagenomic and metatranscriptomic gene copies between pre-drought and drought conditions, PCA was performed using the built-in R ‘prcomp’ function on log-transformed NMR data for all metabolites, FTICR data collapsed at the compound class level and VST-transformed normalized gene copy data using the prcomp function in base R and plotted with ggplot2 (v.3.4.1)103. LME models were created as described above but with only ‘Site’ as a random variable to detect significant differences in all NMR metabolites and FTICR compounds classes between pre-drought and drought conditions. To determine which genes were up- or downregulated during drought, DESeq2 was used to find the log2FC and associated P values (false discovery rate (FDR) corrected) between pre-drought vs drought conditions while controlling for differences between ‘Sites’.
For WGCNA eigengene expression and taxonomic composition, LME models were created with only ‘Site’ as a random variable to determine significant differences between (1) pre-drought vs drought conditions across all timepoints and (2) 6 and 48 vs 0 h within each condition (pre-drought or drought, using ‘Time’ in place of ‘Condition’ in the model). P values reported for Pearson correlations between module eigengenes and 13C-enriched VOC fluxes (acetate-C2 and acetone-C2) were corrected for FDR. Enriched KEGG metabolic pathways within each module were calculated using the ClusterProfiler package (v.3.18.1)104, with the KEGG reference database set to ‘ko’.
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.
The metatranscriptomics and metagenomics sequence data are publicly available through Genbank SRA under the following BioProject IDs: PRJNA980752–PRJNA980834. FTICR, NMR, VOC and CO2 data have been deposited at FigShare (https://doi.org/10.6084/m9.figshare.20334537). Individual raw.xml files for the FTICR data and soil temperature/moisture data are included as Supplementary Data. Source data are provided with this paper.
All code used in data analysis and to create figures is available at www.github.com/linneahkh/SoilPyruvate.
Schimel, J. P. & Schaeffer, S. M. Microbial control over carbon cycling in soil. Front. Microbiol. 3, 00348 (2012).
Insam, H. & Seewald, M. S. A. Volatile organic compounds (VOCs) in soils. Biol. Fertil. Soils 46, 199–213 (2010).
CAS Google Scholar
Penuelas, J. et al. Biogenic volatile emissions from the soil. Plant Cell Environ. 37, 1866–1891 (2014).
CAS PubMed Google Scholar
Honeker, L. K., Graves, K. R., Tfaily, M. M., Krechmer, J. E. & Meredith, L. K. The volatilome: a vital piece of the complete soil metabolome. Front. Environ. Sci. 9, 649905 (2021).
Meredith, L. K. & Tfaily, M. M. Capturing the microbial volatilome: an oft overlooked ‘ome. Trends Microbiol. 30, 622–631 (2022).
CAS PubMed Google Scholar
Koppmann, R. in Hydrocarbons, Oils and Lipids: Diversity, Origin, Chemistry and Fate (ed. Wilkes, H.) 267–277 (Springer, 2010).
Warren, C. R. Response of osmolytes in soil to drying and rewetting. Soil Biol. Biochem. 70, 22–32 (2014).
CAS Google Scholar
Schimel, J., Balser, T. C. & Wallenstein, M. Microbial stress-response physiology and its implications for ecosystem function. Ecology 88, 1386–1394 (2007).
PubMed Google Scholar
Malik, A. A. & Bouskill, N. J. Drought impacts on microbial trait distribution and feedback to soil carbon cycling. Funct. Ecol. 36, 1442–1456 (2022).
CAS Google Scholar
Vardharajula, S. & Sk Z, A. Exopolysaccharide production by drought tolerant Bacillus spp. and effect on soil aggregation under drought stress. J. Microbiol. Biotechnol. Food Sci. 4, 51–57 (2014).
Google Scholar
Sun, S., Lei, H. & Chang, S. X. Drought differentially affects autotrophic and heterotrophic soil respiration rates and their temperature sensitivity. Biol. Fertil. Soils 55, 275–283 (2019).
CAS Google Scholar
Schimel, J. P. Life in dry soils: effects of drought on soil microbial communities and processes. Annu. Rev. Ecol. Evol. Syst. 49, 409–432 (2018).
Google Scholar
Brown, R. W., Chadwick, D. R., Zang, H. & Jones, D. L. Use of metabolomics to quantify changes in soil microbial function in response to fertiliser nitrogen supply and extreme drought. Soil Biol. Biochem. 160, 108351 (2021).
CAS Google Scholar
Bouskill, N. J. et al. Belowground response to drought in a tropical forest soil. I. Changes in microbial functional potential and metabolism. Front. Microbiol. 7, 525 (2016).
PubMed Central PubMed Google Scholar
Bouskill, N. J. et al. Belowground response to drought in a tropical forest soil. II. Change in microbial function impacts carbon composition. Front. Microbiol. 7, 323 (2016).
PubMed Central PubMed Google Scholar
Fang, H. et al. Changes in soil heterotrophic respiration, carbon availability, and microbial function in seven forests along a climate gradient. Ecol. Res. 29, 1077–1086 (2014).
CAS Google Scholar
Li, Y. et al. Effects of biochar application in forest ecosystems on soil properties and greenhouse gas emissions: a review. J. Soils Sediment. 18, 546–563 (2018).
CAS Google Scholar
Su, X. et al. Drought changed soil organic carbon composition and bacterial carbon metabolizing patterns in a subtropical evergreen forest. Sci. Total Environ. 736, 139568 (2020).
CAS PubMed Google Scholar
Peñuelas, J. & Llusià, J. BVOCs: plant defense against climate warming? Trends Plant Sci. 8, 105–109 (2003).
PubMed Google Scholar
Asensio, D., Peñuelas, J., Llusià, J., Ogaya, R. & Filella, I. Interannual and interseasonal soil CO2 efflux and VOC exchange rates in a Mediterranean holm oak forest in response to experimental drought. Soil Biol. Biochem. 39, 2471–2484 (2007).
CAS Google Scholar
Bourtsoukidis, E. et al. Strong sesquiterpene emissions from Amazonian soils. Nat. Commun. 9, 2226 (2018).
CAS PubMed Central PubMed Google Scholar
Jardine, K. et al. Dimethyl sulfide in the Amazon rain forest. Glob. Biogeochem. Cycles 29, 19–32 (2015).
CAS Google Scholar
Reichstein, M. et al. Climate extremes and the carbon cycle. Nature 500, 287–295 (2013).
CAS PubMed Google Scholar
Jentsch, A. & Beierkuhnlein, C. Research frontiers in climate change: effects of extreme meteorological events on ecosystems. C. R. Geosci. 340, 621–628 (2008).
Google Scholar
Bore, E. K., Apostel, C., Halicki, S., Kuzyakov, Y. & Dippold, M. A. Soil microorganisms can overcome respiration inhibition by coupling intra- and extracellular metabolism: 13C metabolic tracing reveals the mechanisms. ISME J. 11, 1423–1433 (2017).
CAS Google Scholar
Dijkstra, P. et al. Modeling soil metabolic processes using isotopologue pairs of position-specific 13C-labeled glucose and pyruvate. Soil Biol. Biochem. 43, 1848–1857 (2011).
CAS Google Scholar
Dijkstra, P. et al. Probing carbon flux patterns through soil microbial metabolic networks using parallel position-specific tracer labeling. Soil Biol. Biochem. 43, 126–132 (2011).
CAS Google Scholar
Albright, M. B. N. et al. Differences in substrate use linked to divergent carbon flow during litter decomposition. FEMS Microbiol. Ecol. 96, fiaa135 (2020).
CAS PubMed Google Scholar
Werner, C. et al. Ecosystem fluxes during drought and recovery in an experimental forest. Science 374, 1514–1518 (2021).
CAS PubMed Google Scholar
Drake, H. L., Küsel, K. & Matthies, C. Ecological consequences of the phylogenetic and physiological diversities of acetogens. Antonie Van Leeuwenhoek 81, 203–213 (2002).
CAS PubMed Google Scholar
Akinyede, R., Taubert, M., Schrumpf, M., Trumbore, S. & Küsel, K. Rates of dark CO2 fixation are driven by microbial biomass in a temperate forest soil. Soil Biol. Biochem. 150, 107950 (2020).
CAS Google Scholar
Miao, W. et al. Protist 10,000 Genomes Project. Innovation 1, 100058 (2020).
PubMed Central PubMed Google Scholar
Chen, Y., Zhang, Y., Siewers, V. & Nielsen, J. Ach1 is involved in shuttling mitochondrial acetyl units for cytosolic C2 provision in Saccharomyces cerevisiae lacking pyruvate decarboxylase. FEMS Yeast Res. 15, fov015 (2015).
PubMed Google Scholar
Buu, L.-M., Chen, Y.-C. & Lee, F.-J. S. Functional characterization and localization of acetyl-coA hydrolase, Ach1p, in Saccharomyces cerevisiae. J. Biol. Chem. 278, 17203–17209 (2003).
Fleck, C. B. & Brock, M. Re-characterisation of Saccharomyces cerevisiae Ach1p: fungal coA-transferases are involved in acetic acid detoxification. Fungal Genet. Biol. 46, 473–485 (2009).
CAS PubMed Google Scholar
Carman, A. J., Vylkova, S. & Lorenz, M. C. Role of acetyl coenzyme A synthesis and breakdown in alternative carbon source utilization in Candida albicans. Eukaryot. Cell 7, 1733–1741 (2008).
CAS PubMed Central PubMed Google Scholar
Williams, F. R., Robert Williams, F. & Hager, L. P. Crystalline flavin pyruvate oxidase from Escherichia coli. Arch. Biochem. Biophys. 116, 168–176 (1966).
CAS PubMed Google Scholar
Schreiner, M. E. & Eikmanns, B. J. Pyruvate:quinone oxidoreductase from Corynebacterium glutamicum: purification and biochemical characterization. J. Bacteriol. 187, 862–871 (2005).
CAS PubMed Central PubMed Google Scholar
Carter, K. & Gennis, R. B. Reconstitution of the ubiquinone-dependent pyruvate oxidase system of Escherichia coli with the cytochrome o terminal oxidase complex. J. Biol. Chem. 260, 10986–10990 (1985).
CAS PubMed Google Scholar
Awang, G. M., Jones, G. A. & Ingledew, W. M. The acetone-butanol-ethanol fermentation. Crit. Rev. Microbiol. 15, S33–S67 (1988).
PubMed Google Scholar
Maddox, I. S. The acetone-butanol-ethanol fermentation: recent progress in technology. Biotechnol. Genet. Eng. Rev. 7, 189–220 (1989).
CAS PubMed Google Scholar
Karl, T., Striednig, M., Graus, M., Hammerle, A. & Wohlfahrt, G. Urban flux measurements reveal a large pool of oxygenated volatile organic compound emissions. Proc. Natl Acad. Sci. USA 115, 1186–1191 (2018).
CAS PubMed Central PubMed Google Scholar
Branen, A. L. & Keenan, T. W. Biosynthesis of α-acetolactate and its conversion to diacetyl and acetoin in cell-free extracts of Lactobacillus casei. Can. J. Microbiol. 18, 479–485 (1972).
Kakumanu, M. L., Ma, L. & Williams, M. A. Drought-induced soil microbial amino acid and polysaccharide change and their implications for C-N cycles in a climate change world. Sci. Rep. 9, 10968 (2019).
PubMed Central PubMed Google Scholar
Kempf, B. & Bremer, E. Uptake and synthesis of compatible solutes as microbial stress responses to high-osmolality environments. Arch. Microbiol. 170, 319–330 (1998).
CAS PubMed Google Scholar
Trchounian, A. A. A direct interaction between the H+-F1F0-ATPase and the K+ transport within the membrane of anaerobically grown bacteria. Bioelectrochem. Bioenerg. 33, 1–10 (1994).
Abdel-Hamid, A. M., Attwood, M. M. & Guest, J. R. Pyruvate oxidase contributes to the aerobic growth efficiency of Escherichia coli. Microbiology 147, 1483–1498 (2001).
CAS PubMed Google Scholar
Wolfe, A. J. The acetate switch. Microbiol. Mol. Biol. Rev. 69, 12–50 (2005).
CAS PubMed Central PubMed Google Scholar
KEGG Pathway Database (Kanehisa Laboratories,); www.genome.jp/kegg/pathway.html (accessed July 15, 2022).
Barnard, R. L., Osborne, C. A. & Firestone, M. K. Responses of soil bacterial and fungal communities to extreme desiccation and rewetting. ISME J. 7, 2229–2241 (2013).
CAS PubMed Central PubMed Google Scholar
Waring, B. G. & Hawkes, C. V. Short-term precipitation exclusion alters microbial responses to soil moisture in a wet tropical forest. Microb. Ecol. 69, 843–854 (2015).
PubMed Google Scholar
Sun, Y. et al. Drought stress induced increase of fungi:bacteria ratio in a poplar plantation. CATENA 193, 104607 (2020).
Börjesson, T., Stöllman, U. & Schnürer, J. Volatile metabolites produced by six fungal species compared with other indicators of fungal growth on cereal grains. Appl. Environ. Microbiol. 58, 2599–2605 (1992).
PubMed Central PubMed Google Scholar
Robinson, P. M. & Garrett, M. K. Identification of volatile sporostatic factors from cultures of Fusarium oxysporum. Trans. Br. Mycol. Soc. 52, 293–299 (1969).
CAS Google Scholar
McBride, S. G., Osburn, E. D., Barrett, J. E. & Strickland, M. S. Volatile methanol and acetone additions increase labile soil carbon and inhibit nitrification. Biogeochemistry 145, 127–140 (2019).
Google Scholar
Slessarev, E. W. & Schimel, J. P. Partitioning sources of CO2 emission after soil wetting using high-resolution observations and minimal models. Soil Biol. Biochem. 143, 107753 (2020).
CAS Google Scholar
Pegoraro, E., Rey, A., Abrell, L., Haren, J. & Lin, G. Drought effect on isoprene production and consumption in Biosphere 2 tropical rainforest. Glob. Change Biol. 12, 456–469 (2006).
Google Scholar
van Haren, J. L. M. et al. Drought-induced nitrous oxide flux dynamics in an enclosed tropical forest. Glob. Change Biol. 11, 1247–1257 (2005).
Google Scholar
Pegoraro, E. et al. Effect of elevated CO2 concentration and vapour pressure deficit on isoprene emission from leaves of Populus deltoides during drought. Funct. Plant Biol. 31, 1137–1147 (2004).
CAS PubMed Google Scholar
Leigh, L. S., Burgess, T., Marino, B. D. V. & Wei, Y. D. Tropical rainforest biome of Biosphere 2: structure, composition and results of the first 2 years of operation. Ecol. Eng. 13, 65–93 (1999).
Google Scholar
Fasbender, L., Yáñez-Serrano, A. M., Kreuzwieser, J., Dubbert, D. & Werner, C. Real-time carbon allocation into biogenic volatile organic compounds (BVOCs) and respiratory carbon dioxide (CO2) traced by PTR-TOF-MS, 13CO2 laser spectroscopy and 13C-pyruvate labelling. PLoS ONE 13, e0204398 (2018).
PubMed Central PubMed Google Scholar
Honeker, L. K. et al. Elucidating drought-tolerance mechanisms in plant roots through 1H NMR metabolomics in parallel with MALDI-MS, and NanoSIMS imaging techniques. Environ. Sci. Technol. 56, 2021–2032 (2022).
CAS PubMed Google Scholar
Ladd, S. N. et al. Leaf-level metabolic changes in response to drought affect daytime CO2 emission and isoprenoid synthesis. Preprint at bioRxiv https://doi.org/10.1101/2022.04.29.490001 (2022).
Deming, B. L. et al. Measurements of delays of gas-phase compounds in a wide variety of tubing materials due to gas–wall interactions. Atmos. Meas. Tech. 12, 3453–3461 (2019).
CAS Google Scholar
Holzinger, R. PTRwid: a new widget tool for processing PTR-TOF-MS data. Atmos. Meas.Tech. https://doi.org/10.5194/amt-8-3903-2015 (2015).
de Gouw, J. et al. Sensitivity and specificity of atmospheric trace gas detection by proton-transfer-reaction mass spectrometry. Int. J. Mass Spectrom. 223–224, 365–382 (2003).
Google Scholar
Holzinger, R. et al. Validity and limitations of simple reaction kinetics to calculate concentrations of organic compounds from ion counts in PTR–MS. Atmos. Meas.Tech. https://doi.org/10.5194/amt-12-6193-2019 (2019).
de Gouw, J. & Warneke, C. Measurements of volatile organic compounds in the earth’s atmosphere using proton-transfer-reaction mass spectrometry. Mass Spectrom. Rev. 26, 223–257 (2007).
PubMed Google Scholar
Krechmer, J. et al. Evaluation of a new vocus reagent-ion source and focusing ion-molecule reactor for use in proton-transfer-reaction mass spectrometry. Anal. Chem. 90, 12011–12018 (2018).
Claflin, M. S. et al. An in situ gas chromatograph with automatic detector switching between PTR- and EI-TOF-MS: isomer-resolved measurements of indoor air. Atmos. Meas. Tech. 14, 133–152 (2021).
Gil-Loaiza, J. et al. Versatile soil gas concentration and isotope monitoring: optimization and integration of novel soil gas probes with online trace gas detection. Biogeosciences 19, 165–185 (2022).
Braun, A. et al. Reviews and syntheses: Heterotrophic fixation of inorganic carbon – significant but invisible flux in environmental carbon cycling. Biogeosciences 18, 3689–3700 (2021).
Pugliese, G. et al. The effect of prolonged drought and recovery on soil VOC fluxes in an experimental rainforest. Nat. Commun. (in the press).
Koebsch, F., Glatzel, S. & Jurasinski, G. Vegetation controls methane emissions in a coastal brackish fen. Wetl. Ecol. Manage. 21, 323–337 (2013).
CAS Google Scholar
Jurasinski, G., Koebsch, F. & Hagemann, U. Flux rate calculation from dynamic closed chamber measurements. R package v0.3-0.1 (2022).
Tfaily, M. M. et al. Advanced solvent based methods for molecular characterization of soil organic matter by high-resolution mass spectrometry. Anal. Chem. 87, 5206–5215 (2015).
CAS PubMed Google Scholar
Tfaily, M. M. et al. Single-throughput complementary high-resolution analytical techniques for characterizing complex natural organic matter mixtures. J. Vis. Exp. https://doi.org/10.3791/59035 (2019).
Tfaily, M. M. et al. Sequential extraction protocol for organic matter from soils and sediments using high resolution mass spectrometry. Anal. Chim. Acta 972, 54–61 (2017).
CAS PubMed Google Scholar
Dittmar, T., Koch, B., Hertkorn, N. & Kattner, G. A simple and efficient method for the solid-phase extraction of dissolved organic matter (SPE-DOM) from seawater. Limnol. Oceanogr. Methods 6, 230–235 (2008).
Willcott, M. R. MestRe Nova. J. Am. Chem. Soc. 131, 13180 (2009).
CAS Google Scholar
Tolić, N. et al. Formularity: software for automated formula assignment of natural and other organic matter from ultrahigh-resolution mass spectra. Anal. Chem. 89, 12659–12665 (2017).
PubMed Google Scholar
Tfaily, M. M. et al. Vertical stratification of peat pore water dissolved organic matter composition in a peat bog in Northern Minnesota. J. Geophys. Res. Biogeosci. 123, 479–494 (2018).
CAS Google Scholar
Ayala-Ortiz, C. O. et al. MetaboDirect: an analytical pipeline for the processing of FTICR-MS-based metabolomics data. Microbiome 11, 38 (2023).
Clum, A. et al. DOE JGI metagenome workflow. mSystems 6, e00804–e00820 (2021).
CAS PubMed Central PubMed Google Scholar
Bushnell, B., Rood, J. & Singer, E. BBMerge-accurate paired shotgun read merging via overlap. PLoS One 12, e0185056 (2017).
Nurk, S., Meleshko, D., Korobeynikov, A. & Pevzner, P. A. metaSPAdes: a new versatile metagenomic assembler. Genome Res. 27, 824–834 (2017).
CAS PubMed Central PubMed Google Scholar
Li, D., Liu, C.-M., Luo, R., Sadakane, K. & Lam, T.-W. MEGAHIT: an ultra-fast single-node solution for large and complex metagenomics assembly via succinct de Bruijn graph. Bioinformatics 31, 1674–1676 (2015).
CAS PubMed Google Scholar
Lowe, T. M. & Eddy, S. R. tRNAscan-SE: a program for improved detection of transfer RNA genes in genomic sequence. Nucleic Acids Res. 25, 955–964 (1997).
CAS PubMed Central PubMed Google Scholar
Nawrocki, E. P. & Eddy, S. R. Infernal 1.1: 100-fold faster RNA homology searches. Bioinformatics 29, 2933–2935 (2013).
CAS PubMed Central PubMed Google Scholar
Bland, C. et al. CRISPR recognition tool (CRT): a tool for automatic detection of clustered regularly interspaced palindromic repeats. BMC Bioinformatics 8, 209 (2007).
PubMed Central PubMed Google Scholar
Hyatt, D. et al. Prodigal: prokaryotic gene recognition and translation initiation site identification. BMC Bioinformatics 11, 119 (2010).
PubMed Central PubMed Google Scholar
Lomsadze, A., Gemayel, K., Tang, S. & Borodovsky, M. Modeling leaderless transcription and atypical genes results in more accurate gene prediction in prokaryotes. Genome Res. 28, 1079–1089 (2018).
Kanehisa, M. & Sato, Y. KEGG Mapper for inferring cellular functions from protein sequences. Protein Sci. 29, 28–35 (2020).
CAS PubMed Google Scholar
Kiełbasa, S. M., Wan, R., Sato, K., Horton, P. & Frith, M. C. Adaptive seeds tame genomic sequence comparison. Genome Res. 21, 487–493 (2011).
PubMed Central PubMed Google Scholar
R Core Team. R: A Language and Environment for Statistical Computing (R Foundation for Statistical Computing) (2020).
Pinheiro, J., Bates, D., DebRoy, S., Sarkar, D. & R Core Team Linear and Nonlinear Mixed Effects Models. R package v.3 (2007).
Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014).
Langfelder, P. & Horvath, S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9, 559 (2008).
PubMed Central PubMed Google Scholar
Langfelder, P. & Horvath, S. Fast R functions for robust correlations and hierarchical clustering. J. Stat. Softw. 46, i11 (2012).
PubMed Central PubMed Google Scholar
Mevik, B.-H. & Wehrens, R. The pls package: principal component and partial least squares regression in R. J. Stat. Softw. https://doi.org/10.18637/jss.v018.i02 (2007).
Lu, B., Castillo, I., Chiang, L. & Edgar, T. F. Industrial PLS model variable selection using moving window variable importance in projection. Chemometr. Intell. Lab. Syst. 135, 90–109 (2014).
CAS Google Scholar
Liland, K. H., Mehmood, T. & Sæbø, S. plsVarSel: Variable Selection in Partial Least Squares. R package v.0.9.6 (2020).
Wilkinson, L. ggplot2: elegant graphics for data analysis by WICKHAM, H. Biometrics 67, 678–679 (2011).
Google Scholar
Wu, T. et al. clusterProfiler 4.0: a universal enrichment tool for interpreting omics data. Innovation 2, 100141 (2021).
CAS PubMed Central PubMed Google Scholar
Download references
A portion of this research was performed under the Facilities Integrating Collaborations for User Science (FICUS) programme (proposal: https://doi.org/10.46936/fics.proj.2019.50971/60000130) awarded to M.M.T. and used resources at the DOE Joint Genome Institute (https://ror.org/04xm1d337) and the Environmental Molecular Sciences Laboratory (https://ror.org/04rc0xn13), which are DOE Office of Science User Facilities operated under Contract Nos. DE-AC02-05CH11231 (JGI) and DE-AC05-76RL01830 (EMSL). This material was based upon work supported by the US Department of Energy, Office of Science, Biological and Environmental Research Program under Award Number DE-SC0023189 awarded to L.K.M. This research was also supported in part by the European Research Council (ERC; Grant Number 647008) awarded to C.W., Department of Energy, Office of Science Biological and Environmental Research Grant (DE-SC0021349), and NSF CAREER Award (No. 2045332) to L.K.M. L.K.H. was supported by Biosphere 2 through the office of the Senior Vice President for Research Innovation and Impact at the University of Arizona, NSF CAREER Award (No. 2045332), and the BIO5 Postdoctoral Fellowship. G.P. was supported by the German Federal Ministry of Education and Research (BMBF contract 01LB1001A – ATTO+) and the Max Planck Society. We Also acknowledge financial support from the Philecology Foundation.
Jane Fudyma
Present address: Department of Plant Pathology, University of California, Davis, CA, USA
Jordan E. Krechmer
Present address: Bruker Daltonics Inc., Billerica, MA, USA
Eva Y. Pfannerstill
Present address: Department of Environmental Science, Policy, and Management, University of California, Berkeley, CA, USA
Biosphere 2, University of Arizona, Tucson, AZ, USA
Linnea K. Honeker & Laura K. Meredith
School of Natural Resources and the Environment, University of Arizona, Tucson, AZ, USA
Linnea K. Honeker, Juliana Gil-Loaiza & Laura K. Meredith
Ecosystem Physiology, Faculty of Environment and Natural Resources, University of Freiburg, Freiburg, Germany
Giovanni Pugliese, Johannes Ingrisch, L. Erik Daber, Jürgen Kreuzwieser, S. Nemiah Ladd & Christiane Werner
Max Planck Institute for Chemistry, Atmospheric Chemistry Department, Mainz, Germany
Giovanni Pugliese, Eva Y. Pfannerstill & Jonathan Williams
Department of Ecology, Universität Innsbruck, Innsbruck, Austria
Johannes Ingrisch & Kathiravan Meeran
Department of Environmental Sciences, University of Arizona, Tucson, AZ, USA
Jane Fudyma, Gina Hildebrand, Christian Ayala-Ortiz, Viviana Freire-Zapata & Malak M. Tfaily
Joint Genome Institute, Walnut Creek, CA, USA
Elizabeth Carpenter & Esther Singer
Geo-Biosphere Interactions, Department of Geosciences, University of Tuebingen, Tuebingen, Germany
Lingling Shi & Michaela A. Dippold
Environmental Molecular Science Laboratory (EMSL), Earth and Biological Sciences Division, Pacific Northwest National Laboratory, Richland, WA, USA
David W. Hoyt, Rosalie K. Chu & Jason Toyoda
Aerodyne Research, Inc., Billerica, MA, USA
Jordan E. Krechmer & Megan S. Claflin
Department of Environmental Sciences, University of Basel, Basel, Switzerland
S. Nemiah Ladd
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
You can also search for this author in PubMed Google Scholar
C.W., L.K.M. and S.N.L. conceived the project and experimental design. J.I., J.F., J.G.-L., J.E.K., M.A.D., L.E.D., S.N.L., C.W., M.M.T. and L.K.M. contributed to experimental design. J.I., J.F., J.G.-L., J.E.K., J.K., M.S.C., L.E.D., M.A.D., E.Y.P., L.S., S.N.L and L.K.M. performed the experiments in the field. J.E.K., J.W., J.K., D.W.H., R.K.C., J.T. and M.M.T. contributed equipment and analysis tools. E.C., E.S., L.K.H., G.H. and D.W.H. contributed to sample preparation and processing. L.K.H., G.P., J.E.K., K.M., J.I., D.W.H., C.A.-O. and V.F.-Z. analysed the data. L.K.H., G.P. and J.E.K. contributed to making figures. L.K.H., L.K.M, S.N.L. and M.M.T. wrote the paper.
Correspondence to Laura K. Meredith.
The authors declare no competing interests.
Nature Microbiology thanks Riikka Rinnan, Steve McBride, Steffen Kolb and the other, anonymous, reviewer(s) for their contribution to the peer review of this work.
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
a) Topological map of Biosphere 2 Tropical Rainforest, with locations of Sites 1 − 3. All circles represent collar locations where CO2 and VOC data were collected and large circles represent a sub selection of sites collars where additional soil samples were collected for metatranscriptomics, metagenomics, and metabolomics. b) Diagram showing soil pyruvate experimental setup with stencil placed inside a soil collar where 13C1- or 13C2-pyruvate was added within each 1 cm×1 cm square. An automatic soil chamber was placed over the soil collar during the experiment and continuous measurements of 12C and 13C-CO2 and VOCs were measured.
a) Total flux of CO2 during pre-drought and drought conditions. Each point represents a single measurement per soil chamber (~ 64 measurements), of which there were 3 replicates each that were injected with either 13C1- or 13C2-pyruvate for each site (total n = 18 each for pre-drought and drought), and lines show the data smoothed with the surrounding shaded area showing ± SEM. B) Proportion of C from pyruvate allocated to biosynthesis. This was calculated as 13C-CO2-C1/ (13C-CO2-C1 + 13C-CO2-C2) on each set of C1/C2 chambers per site (n = 9 each for pre-drought and drought) with continuous emission data binned to 3 or 6 h intervals, from 0 to 48 h post pyruvate injection. (n = 9 each for pre-drought and drought). Boxes represent Q1 - Q3 with center line indicating median, and bars extending to maximum and minimum values, excluding outliers.
Source data
Two nearby locations a) and b) had two peaks with retention times 217 s and 275 s for C4H6O2 corresponding to two different compounds. The retention time of 217 s matched the expected retention time for diacetyl, but the retention time of 275 s did not match with any known compounds. Therefore, C4H6O2 may represent diacetyl and/or one additional unidentifiable compound, and is referred to as diacetyl+ in the manuscript. The diacetyl peak (retention time of 217 s) is indicated with an arrow.
PCA of (a) metatranscriptomic and (b) metagenomics data. Taxonomic profiles as inferred from metatranscriptomics data at the (c) kingdom-level for active archaea, bacteria, eukaryota, and viruses, d) phylum-level for bacteria and archaea, and e) phylum-level for fungi. Taxonomic profiles as inferred from metagenomics data at the f) kingdom-level for total archaea, bacteria, eukaryota, and viruses, g) phylum-level for bacteria and archaea, and h) phylum-level for fungi. For taxonomic profiles, arrows indicate direction of changes in relative abundance (up for increased, and down for decreased relative abundance in drought compared to pre-drought conditions). Arrows with no stars means the phylum was present in one condition (pre-drought or drought) and absent in the other condition. For panels a-b, P-values were determined using PERMANOVA on Bray-Curtis distance matrices. For panels c-h, P-values were determined using Linear Mixed Effect Models. *, P < 0.05; **, P < 0.01; (exact P-values: P = 0.010 [Proteobacteria], P = 0.016 [Gemmatimonadetes], P = 0.0012 [Actinobacteria] (panel b), P = 0.0050 [Ascomycota] (panel c)). metaT, metatranscriptomics; metaG, metagenomics; B, bacteria; A, archaea; F, fungi.
Source data
Log2-fold change (FC) of genes involved in cycling of acetate, acetone, and diacetyl+ at 6 and 48 h in relation to 0 h to show effect of pyruvate addition on microbial activity during a) pre-drought and b) drought conditions. Log2-FC values to indicate up-or down-regulation calculated with DESeq2. Production vs. consumption indicated by black and grey bars, respectively. Significantly up- or down-regulated genes are outlined in color (blue = upregulated, red = downregulated), with exact P-values from DESeq2 analysis shown.); *, P < 0.05; **, P < 0.01; ***, P < 0.001.
Source data
Gene copies (transcripts) of a) acetate- and b) acetone-cycling genes at 0, 6, and 48 h post pyruvate injection. Bars represent Q1 - Q3 with center line indicating median, and bars extending to maximum and minimum values, excluding outliers, across all chambers and sites (n = 9 each for pre-drought and drought).
Source data
a) PCA biplot of samples collected during pre-drought (green; 0 h [n = 5]) or drought (orange; 0 h [n = 5] and 6 h [n = 1]) showing relationships between metabolic profiles identified with NMR which were mostly identified as primary metabolites with relatively low molecular weight. Arrows represent the loadings of individual metabolites. b) PCA biplot of samples collected during pre-drought (green) or drought (orange) showing relationships between profiles of metabolic classes identified with FTICR-MS, capturing mostly high molecular weight metabolites including secondary metabolites. Arrows represent the loadings of metabolic classes driving these patterns. P-values indicate significance of clustering between pre-drought and drought conditions using PERMANOVA on bray-curtis dissimilarity matrix. c) Each stacked bar represents the composition of compound classes identified with FTICR-MS for all samples (n = 18 [6 chambers sampled at 0, 6, and 48 h post pyruvate addition] averaged across pre-drought and drought conditions). Exact P-values are shown as determined by Linear Mixed Effects Models; *, P < 0.05; **, P < 0.01; ***, P < 0.001.
Source data
a) WGCNA cluster dendrogram reveals 9 clusters of co-expressed genes, or modules. The four modules with significant correlations to drought or pre-drought (pink, green, magenta, and brown) are indicated. b) Correlations of the pink, green, magenta, and brown module eigengenes with environmental conditions (pre-drought = 0, drought = 1) and efflux of 13C-enriched (VOC 13C/(12C + 13C)) acetate and acetone from chambers receiving 13C2-pyruvate (acetate-C2 and acetone-C2, respectively) across three time points (0, 6, and 48 h post pyruvate injection; total n = 27 [3 measurements each for each soil collar receiving 13C2-pyruvate]). Pearson correlation coefficients are shown with FDR-corrected p-values below in parentheses. White boxes indicate no correlation (P < 0.1). *, P < 0.05; **, P < 0.01; ***, P < 0.001.
Source data
Supplementary Figs.1–8, and Tables 1 and 2.
Folder containing all raw XML data for FTICR.
Soil moisture and temperature data.
dCO2 flux, acetate, acetone, diacetyl flux and area under the curve data, and CO2 C1-C2 differences.
DESeq2 results of VOC cycling genes in Fig. 3.
WGCNA eigengene expression for pink, red, green and magenta modules.
CO2 flux and CO2 C1-C2 ratios.
metaT/metaG PCA coordinates and metaT/metaG taxonomy.
DESeq2 results of VOC cycling genes in Extended Data Fig. 5.
Gene/transcript count table, sample metadata and trait data for WGCNA.
NMR data and FTICR class compounds.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
Reprints and Permissions
Honeker, L.K., Pugliese, G., Ingrisch, J. et al. Drought re-routes soil microbial carbon metabolism towards emission of volatile metabolites in an artificial tropical rainforest. Nat Microbiol 8, 1480–1494 (2023). https://doi.org/10.1038/s41564-023-01432-9
Download citation
Received: 08 July 2022
Accepted: 19 June 2023
Published: 31 July 2023
Issue Date: August 2023
DOI: https://doi.org/10.1038/s41564-023-01432-9
Anyone you share the following link with will be able to read this content:
Sorry, a shareable link is not currently available for this article.
Provided by the Springer Nature SharedIt content-sharing initiative