Bacterial community assembly on tobacco leaves was affected by fermentation time and leaf position

After quality filtration, a total of 11,011,441 high-quality sequences were retained in 200 samples and these reads were assigned into 10,227 ASVs. The ASV table was then normalized to 24,000 reads in each sample for the analysis of microbial diversity and composition. PCoA ordinations and PERMANOVA analysis showed that the assembly of bacterial community on tobacco leaves during the fermentation process was determined by fermentation time, leaf position and fermentation site to varying degrees (Fig. 1). The PCoA plot showed the two axes (PCo1 and PCo2) explained 42.82% of the variance and the significant dissimilarity of the whole bacterial community (Figs. 1A and S1). The results of PERMANOVA revealed that fermentation time explained the largest variation in the whole bacterial community (49.54%), followed by leaf position (7.24%) and fermentation site (2.25%) (Fig. 1A and Table S1). In the early stage of fermentation (M0, M3 and M6), the assembly of bacterial community was mainly explained by leaf position (34.30%, 43.36% and 66.07%, respectively). The relative contribution of fermentation site in shaping microbial community was varied (13.41–41.60%) with fermentation stages (Fig. 1B and Table S2).

Fig. 1

The changes of bacterial community diversity on tobacco leaves during the fermentation process. A Principal coordinate analysis (PCoA) showing the effects of fermentation time, leaf position and fermentation site on the structures of bacterial community on tobacco leaves. B PCoA showing the effects of leaf position and fermentation site on the structures of bacterial community at different fermentation stages. Asterisks denote significant differences between samples from top and middle leaves (one-way ANOVA followed by a post hoc Tukey test, ***P < 0.001). C Shannon index and Chao1 showing the effects of fermentation time and leaf position on the α-diversity of bacterial community on tobacco leaves

Four alpha diversity indices including Shannon index, Chao1 index, the observed OTUs and Simpson diversity were calculated (Figs. 1C, S2 and S3). It was found that the bacterial alpha diversity was significantly affected (Kruskal-Wallis; P < 0.01) by fermentation time. With the progress of fermentation, the bacterial alpha diversity showed a trend of reducing first (M3, M6 and M9) and then increasing (M12) (Fig. S3). For different leaf positions, we found the community diversity of top leaves was significantly higher than that of middle leaves at the beginning of fermentation. However, the difference was gradually reduced or even disappeared in the fermentation process (Figs. 1C and S2).

The calculation results of deterministic (|βNTI|≥ 2) and stochastic (|βNTI|< 2) processes showed that the βNTI values decreased significantly during the first 9 months of the fermentation process (P < 0.01) (Fig. 2A), indicating that bacterial community assembly altered from stochastic to deterministic processes. Null model analysis was further performed to seek the community assembly processes over different fermentation stages (Fig. 2B). A higher relative contribution of stochastic processes primarily associated with homogenizing dispersal and undominated (drift and/or diversification) was found in M0 (91.73% for middle leaves and 81.93% for top leaves) and M12 (80.99% for middle leaves and 86.51% for top leaves) bacterial communities. Conversely, bacterial communities in M3 (88.22% for middle leaves and 85.56% for top leaves), M6 (88.45% for middle leaves and 92.46% for top leaves) and M9 (86.91% for middle leaves and 89.56% for top leaves) were driven by deterministic processes mainly belonging to homogeneous selection. Totally, stochastic process was the primary contributor to bacterial community assembly at the early (M0) and late (M12) stages and deterministic process exerted a greater influence on community in the middle stages (M3, M6 and M9) of tobacco fermentation.

Fig. 2

Deterministic and stochastic processes in bacterial community assembly during the fermentation process. A Relative contribution of determinism and stochasticity on bacterial community assembly on tobacco leaves during the fermentation process based on the β-Nearest Taxon Index (βNTI) values. B The relative importance of five ecological processes on bacterial community at different fermentation stages

Different marker microbes were detected at different fermentation stages

Taxonomic classification showed that the bacterial community on fermented tobacco leaves was mainly composed of Proteobacteria (relative abundance ranging from 39.27 to 94.97%) and Firmicutes (relative abundance ranging from 2.86 to 59.10%) at the phylum level (Fig. 3A). In the early stages of fermentation (M0, M3 and M6), the content of Firmicutes was more abundant on middle leaves than that on top leaves, but Proteobacteria presented a contrary trend (Fig. 3A). The taxonomic composition at the genus level showed a significant variation with the fermentation process (Fig. 3B). LefSe analysis was conducted at the genus level to determine the enrichment/depletion pattern between middle and top leaves (Fig. S4) or among different fermentation time (Fig. 4) with an LDA score threshold of 4. We found that the relative abundance of Terribacillus and Bacillus on middle leaves were much higher in most fermentation stages (M0, M3, M6 and M9), but Sphingomonas and Pantoea were more abundant on top leaves (Fig. S4), which indicated the potential marker genera on middle and top tobacco leaves during fermentation.

Fig. 3

The changes of bacterial composition on tobacco leaves during the fermentation process. A The changes of bacterial composition at the phylum level. B The changes of bacterial composition at the genus level

Considering the constant changes of microorganisms during fermentation, the specific microorganisms of middle and top leaves in different fermentation stages were excavated at the genus level. Compared with other fermentation stages, the genera enriched on middle (Enterobacter and Kosakonia) and top leaves (Sphingomonas, Stenotrophomons, Methylobacterium-Methylorubrum and Pseudocitrobacter) in M0 were significantly different (Fig. 4A). However, as the fermentation process progressed, the microorganisms enriched in the middle and top leaves tended to be consistent in the same fermentation stage, indicating that the enriched microorganisms may play an important role in the corresponding stage. For example, Franconibacter was enriched in M6, Pseudomonas was enriched in M9 and Xylella was enriched in M12, for both middle and top leaves. Besides, Terribacillus and Pantoea were enriched in the early stage of fermentation, while the enrichment periods in middle and top leaves were different (Fig. 4A).

Differential analysis of ASVs confirmed the constantly variation of microorganisms (Fig. 4B). Compared with the original bacterial community (M0), there were 63, 67, 72 and 88 ASVs enriched in M3, M6, M9 and M12 for middle leaves, and 48, 60, 49 and 85 ASVs enriched in M3, M6, M9 and M12 for top leaves. Notably, a large number of ASVs were depleted in M3 and M6 for both middle (143 and 128, respectively) and top (215 and 190, respectively) leaves when compared with those in M0, which indicated the selection and filtration for microorganisms on the surface of tobacco leaves during the fermentation process.

Fig. 4

The enrichment and depletion patterns of bacterial composition during the fermentation process. A The LEfSe analysis with an LDA score threshold of 4 showing the specific bacterial genera of middle and top tobacco leaves at different fermentation stages. B The volcano plots showing the enriched or depleted bacterial ASVs of middle and top tobacco leaves at different fermentation stages (M3, M6, M9 and M12) compared with the original bacterial community (M0)

The complexity and stability of microbial co-occurrence networks were changed during fermentation

To assess the impact of fermentation time on potential interaction dynamics within the bacterial community, co-occurrence networks between bacteria at different fermentation stages were constructed (Fig. 5). The results showed that bacterial network patterns shifted clearly across five fermentation stages (Fig. 5A). The complexity (nodes and edges) and stability (the proportion of negative connections) of bacterial networks of middle and top leaves presented the same trend during the fermentation process, but the networks of top leaves were more connected and complex than those of middle leaves at the same fermentation stage (Fig. 5B). Specifically, bacterial networks of M3, M6 and M9 showed significantly lower connections (edges) than the networks of M0 and M12 for both top and middle leaves (Fig. 5B, Tables S3 and S4). However, the proportion of negative interactions which represented the network stability were significantly higher in M3 (49.21% and 44.99% for middle and top leaves) and M6 (45.04% and 44.21% for middle and top leaves) (Fig. 5B, Tables S3 and S4).

According to the within-module connectivity (Zi) and among-module connectivity (Pi) of nodes, they could be further defined as peripherals, connectors, module hubs, and network hubs (Fig. 5C). Two module hubs (M12) and 14 connectors (2 in M0, 9 in M3 and 3 in M12) were detected in middle leaves, but only two module hubs (1 in M0 and 1 in M9) and two connectors (1 in M0 and 1 in M9) were detected in top leaves (Fig. 5C). All these module hubs and connectors belonged to Alphaproteobacteria (5) or Gammaproteobacteria (14) except for ASV1 (belonging to Bacilli) in M0 (Table S5). The degree of module hubs and connectors in middle leaves ranged from 2 to 25, but that of top leaves was significantly higher (ranging from 12 to 59) (Fig. 5D and E). For middle leaves, we found that 6 connectors (ASV82, ASV50, ASV94, ASV 12, ASV9 and ASV60) belonging to Enterobacteriaceae and Moraxellaceae in M3 and 2 connectors belonging to Sphingomonas in M12 (ASV541 and ASV 86) (Fig. 5D) were enriched members in M3 and M12 (Fig. 4B), respectively. For top leaves, a connector (ASV9) belonging to Enterobacteriaceae (Fig. 5E) was also an enriched member in M9 (Fig. 4B).

Fig. 5

The co-occurrence networks of bacterial community on tobacco leaves during the fermentation process. A Bacterial co-occurrence networks on tobacco leaves during the fermentation process. Nodes represent ASVs. The size of each node represents connections with other nodes and the colors of nodes represent different phyla. The edges between the nodes indicate strong and significant (P < 0.01) correlations. A blue line indicates a positive interaction, while a red line indicates a negative interaction. B Number of positive and negative interactions in bacterial networks at different fermentation stages. C The scatter plot of within-module connectivity (Zi) and among-module connectivity (Pi) showing the topological role of each ASV. D The number and types (positive or negative) of edges showing the interactions of module hubs and connectors on middle leaves with other ASVs. The nodes marked in red were bacteria that are significantly enriched at the corresponding fermentation stage compared with the original bacterial community (M0). E The number and types (positive or negative) of edges showing the interactions of module hubs and connectors on top leaves with other ASVs

The functional metabolisms of microbial community were changed during fermentation

Metagenomes of bacterial communities were predicted using PICRUSt2 and then annotated by referring to the KEGG database. Shannon index based on KOs (KEGG Orthologs) showed that the KO diversity was significantly affected (Kruskal-Wallis; P < 0.01) by fermentation time (Figs. S5 and 6A). The KO diversities of M3, M6 and M9 were significantly lower than that of M0 (Fig. S5). For different leaf positions, we found the KO diversity of M3 and M12 was significantly lower (P < 0.05) in top leaves than that in middle leaves, but KO diversity of M9 showed a contrary trend (Fig. 6A). PCoA analysis for KO category indicated that there were significant differences in community functions at different fermentation stages (Figs. 6B and S6). The results of PERMANOVA revealed that fermentation time (51.68%) and leaf position (11.27%) were important factors in affecting community functions (Fig. 6B).

Fig. 6

PICRUSt predicted metagenome functions of bacterial community on tobacco leaves at KO level. A Functional diversity of KO genes in bacterial community on top and middle leaves at different fermentation stages. Asterisks denote significant differences between samples from top and middle leaves (one-way ANOVA followed by a post hoc Tukey test, *P < 0.05; **P < 0.01). B PCoA analysis showing the effects of fermentation time, leaf position and fermentation site on the distribution of KO functional genes in bacterial community on tobacco leaves. C Heatmap showing the relative abundance of functional genes involved in carbohydrate metabolism, amino acid metabolism, cell motility, nucleotide metabolism and cell community-prokaryotes

Several metabolic pathways including carbohydrate metabolism, amino acid metabolism, cell motility and cell community-prokaryotes were enriched in M3, M6 and M9 (Fig. 6C). However, the specific pathways enriched in M3 and M6 were significantly different from those in M9 (Figs. S7, S8, S9 and S10). For example, for carbohydrate metabolism, genes involved in galactose metabolism (K02080, K10984, K10985, K10986 and K12112), C5-branched dibasic acid metabolism (K18292, K18290 and K18288) and pentose and glucuronate interconversions (K08092, K16849 and K13874) were more abundant in M3 and M6, but genes involved in butanoate metabolism (K01799) and glyoxylate and dicarboxylate metabolism (K00281, K00124 and K00127) were more abundant in M9 (Figs. 6C and S7). For amino caid metabolism, genes involved in arginine and proline metabolism, tryptophan metabolism, histidine metabolism and lysine biosynthesis were more abundant in M9 (Figs. 6C and S8). For cell motility and cell community-prokaryotes, the genes involved bacterial chemotaxis, flagellar assembly, quorum sensing and biofilm formation in M3, M6 and M9 were much higher than M0 and M12 (Figs. 6C, S9 and S10).