Bacterial genome size and gene functional diversity negatively correlate with taxonomic diversity along a pH gradient

Abstract

Bacterial gene repertoires reflect adaptive strategies, contribute to ecosystem functioning and are limited by genome size. However, gene functional diversity does not necessarily correlate with taxonomic diversity because average genome size may vary by community. Here, we analyse gene functional diversity (by shotgun metagenomics) and taxonomic diversity (by 16S rRNA gene amplicon sequencing) to investigate soil bacterial communities along a natural pH gradient in 12 tropical, subtropical, and temperate forests. We find that bacterial average genome size and gene functional diversity decrease, whereas taxonomic diversity increases, as soil pH rises from acid to neutral; as a result, bacterial taxonomic and functional diversity are negatively correlated. The gene repertoire of acid-adapted oligotrophs is enriched in functions of signal transduction, cell motility, secretion system, and degradation of complex compounds, while that of neutral pH-adapted copiotrophs is enriched in functions of energy metabolism and membrane transport. Our results indicate that a mismatch between taxonomic and functional diversity can arise when environmental factors (such as pH) select for adaptive strategies that affect genome size distributions.

Fig. 1. Bacterial community structure along a pH gradient. a Distribution of 12 forests.

Location of 12 forests in China along a latitudinal gradient. GH Genhe, LS Liangshui, CBS Changbaishan, DLS Donglingshan, BTM Baotianman, TTS Tiantongshan, BDGS Badagongshan, GTS Gutianshan, HSD Heishiding, DHS Dinghushan, NG Nonggang, XSBN Xishuangbanna. Color codes for sites in (a–c, e) is consistent with other figures of this paper. The data for map was download from DATAV.GeoAltas (http://datav.aliyun.com/portal/school/atlas/area_selector) and visualized by ggplot2 package (https://ggplot2.tidyverse.org/). b, c Bacterial community composition in association with environmental variables. Principal coordinate (PCo) analysis with environmental fitting (envfit) showing association of 16S rRNA gene amplicon-based bacterial taxonomic composition (b) and metagenome-based bacterial functional composition (c) with soil pH as well as other biotic and abiotic variables (the arrowed lines). The strength and significance of association between PCo vectors and variables are provided in Supplementary Fig. 1. MAP mean annual precipitation, MAT mean annual temperature, TC total carbon, TN total nitrogen, TP total phosphorus, ACa available calcium, AMg available magnesium, AFe available iron, AK available potassium, C_N carbon nitrogen ratio, C_P carbon phosphorus ratio, N_P nitrogen phosphorus ratio. d Bacterial phylum composition along a pH gradient. e Regression curve of detected bacterial phyla against soil pH. Linear regression model with two-sided test was used for the statistical analysis, and adjusted R-squared was used. Both the 16S rRNA gene amplicon (R = −0.674, P = 3.925e-06) and metagenome (R = −0.777, P = 1.52e-08) showed the loss of Acidobacteria with increasing soil pH along the horizontal axis. The relative abundance of Actinobacteria, Planctomycetes, Chloroflexi tends to increase with increasing pH. The relative abundances of Verrucomicrobia, Bacteroidetes, Nitrospirae, Elusimicrobia and Fibrobacteres peak at around pH 5.5. Note that estimates of 16S rRNA gene amplicon and metagenome sequencing are not necessarily consistent. The relative abundance of Proteobacteria remains almost unchanged across the pH gradient. n = 36 samples. The grey area around the smooth line indicates the 95% confidence interval. Source data are provided as a Source Data file.

Fig. 2. Metagenomic traits along a pH gradient.

a Bacterial average genome size (AGS) decreased as pH changed from acid to neutral.Bacterial AGS are detected by analyzing shotgun metagenome using MicrobeCensus pipeline. b Bacterial average 16S rRNA gene copy number (ACN) not significantly associated with soil pH. Bacterial ACN are detected by analyzing shotgun metagenome using the method from ref. . c Bacterial GC content (GC%) increased as pH changed from acid to neutral. Bacterial GC% are detected by analyzing shotgun metagenome using Quast software. d Bacterial estimated growth rate was unaffected by soil pH. Bacterial growth rate (minimal doubling time) was detected by analyzing shotgun metagenome using the gRodon pipeline. Linear regression model with two-sided test was used for the statistical analysis, and adjusted R-squared was used. n = 36 samples. The grey area around the smooth line indicates the 95% confidence interval. Source data are provided as a Source Data file.

Fig. 3. Contrasting distribution patterns for diversities of bacterial taxonomy and functional genes along a pH gradient.

Diversities shown are measured by richness (S), and diversities measured by Shannon’s index (H’) are provided in Supplementary Fig. 7. a Bacterial taxonomic diversity (S.16S) increased as soil pH changed from acid to neutral. Bacterial operational taxonomic units (OTUs) are detected by 16S rRNA gene amplicon metabarcoding sequencing. b Bacterial functional diversity (S.KO) decreased as pH changed from acid to neutral. Bacterial functions are determined from the shotgun metagenome, as annotated by Kyoto Encyclopedia of Genes and Genomes (KEGG) Ontology (KO). c Bacterial diversity of antibiotic resistance genes (S.ARG) decreased as pH changed from acid to neutral. Bacterial antibiotic resistance genes are detected based on shotgun metagenome annotated by the Resfam database. d Bacterial diversity of carbohydrate-active enzymes (S.Cazy) genes decreased as pH changed from acid to neutral. Bacterial carbohydrate-active enzymes genes are detected based on shotgun metagenome annotated by database of CAZy. e Bacterial taxonomic richness (S.16S) negatively correlated with functional gene diversity (S.KO). f–h Bacterial average genome size (AGS) positively correlated with functional diversities as measured by f S.KO, g S.ARG and h S.Cazy. Linear regression model with two-sided test was used for the statistical analysis, and adjusted R-squared was used. n = 36 samples. The grey area around the smooth line indicates the 95% confidence interval. Source data are provided as a Source Data file.

Fig. 4. Functional enrichment and depletion along the pH gradient.

a Network analysis of KOs detected two dominant modules, Neutral pH Module, M1 (blue), containing 2777 vertices, and Acid pH module, M2 (red), containing 4309 vertices. b, c M1 correlates positively and M2 correlates negatively with soil pH. Linear mixed-effects models with two-sided test were used for the statistical analysis. The grey lines are regression lines for each KO, and the colored line are regression line for the average of all KOs. b Neutral pH module (M1, blue), showing that soil pH is significantly, positively correlated with M1 and c Acid pH module (M2, red), showing that soil pH is significantly, negatively correlated with M2. d, e Enrichments of genes in KEGG pathways for Neutral pH (M1, Blue) and Acid pH (M2, Red) modules. Differential expression analysis with two-sided test were used for the statistical analysis. M1 functions were enriched for energy metabolism, membrane transport, citrate cycle, and glyoxylate, dicarboxylate and amino acid metabolism. M2 functions were enriched for bacterial secretion system, cell motility, xenobiotics biodegradation and metabolism, two component systems, metabolism of terpenoids and polyketides, glycan biosynthesis and metabolism, starch and sucrose metabolism, porphyrin metabolism, siderophore and lipopolysaccharide synthesis. f Association between environmental variables and specific genes involved in essential functions drawn from functional pathways (n = 36). In general, soil pH positively correlates with genes in the functions enriched in M1, and negatively correlated with that in M2. Note that associations between environmental variables and genes in specific functions are detailed in Supplementary Figs. 12–18. Source data are provided as a Source Data file.

Fig. 5. Genes involving biogeochemical cycle associated with environmental variables.

a Heatmap showing correlations between environmental variables and genes involving biogeochemical cycling of carbon, nitrogen, phosphorus, sulfur and iron (n = 36 samples). b Seventeen carbon cycle genes positively correlate with soil pH. c Nitrogen cycle genes are positively or negatively correlated with soil pH. d Eight phosphorus cycle genes positively correlate with soil pH, as compared to that two phosphorus cycle genes negatively correlate with soil pH. The diagrams for sulfur and iron cycles can be found in Supplementary Fig. 18. Source data are provided as a Source Data file.

Fig. 6. Conceptual model on the adaptive strategy of soil microbiome along a latitudinal pH gradient across 12 forests.

Oligotrophs characterized by larger genome are adapted to the acidic, resource-poor soil, with an enrichment on functions of cell motility, bacterial chemotaxis, secretion system, signal transduction and complex matter degradation. The copiotrophs characterized by smaller genome are adapted to the neutral, resource-rich soil, with an enrichment on the functions of energy metabolism, membrane transport and amino acid metabolism. Note in our study system the resource availability was poor in acid pH soil, and rich in neutral pH soil. However, resources availability is known to be low in some neutral pH soil and high in some acid pH soil and microbial functional traits and adaptive strategies in these two conditions remain unresolved.