Soil physicochemical parameters and tuber size of T. hemsleyanum
The effects of soil type (stony and non-stony) on the physicochemical parameters of the soil and on the tuber size of T. hemsleyanum were significant (p < 0.05) after two years of growth (Table 1 and Fig. S1 in Additional file 1). Several physicochemical parameters of stony soils, such as electrical conductivity (EC), total carbon (TC), organic matter (OM), and AP were significantly higher (p < 0.05) than those of non-stony soils, regardless of whether bulk soils (BS) or root-zone soils (RZS) were studied.
Table 1
Soil physiochemical properties of the fields associated with non-stony and stony soils. EC, electrical conductivity; TN: total nitrogen; TC, total carbon; OM, organic matter; NH3-N, ammonia nitrogen; NO2-N: nitrite; AP, available phosphate; AK, available potassium; AS, available sulphur. Data are mean ± SD (n = 4 in BS and n = 5 in RZS samples). Different letters indicate a significant difference at p < 0.05.
Samples | Non-BS | St-BS | Non-RZS | St-RZS |
pH | 4.43 ± 0.03b | 4.59 ± 0.07ab | 4.41 ± 0.17b | 4.62 ± 0.15a |
EC µs cm− 1 | 50.6 ± 0.92b | 67.8 ± 6.28a | 50.1 ± 4.34b | 66.8 ± 14.2a |
TN | g kg− 1 | 2.14 ± 0.71a | 2.56 ± 0.04a | 1.95 ± 1.00a | 2.53 ± 0.15a |
TC | 21.6 ± 2.89b | 37.2 ± 2.38a | 23.1 ± 5.13b | 36.9 ± 3.41a |
OM | 33.1 ± 7.28b | 55.3 ± 4.65a | 37.0 ± 12.9b | 55.4 ± 3.58a |
NH3-N | mg kg− 1 | 1.67 ± 0.71a | 1.68 ± 1.02a | 1.18 ± 0.74a | 2.40 ± 1.48a |
NO2-N | 0.04 ± 0.02a | 0.02 ± 0.01a | 0.06 ± 0.05a | 0.013 ± 0.01a |
AP | 0.88 ± 0.17b | 3.08 ± 0.22a | 0.91 ± 0.37b | 2.97 ± 0.91a |
AK | 101.4 ± 5.67a | 115.0 ± 10.2a | 101.9 ± 15.9a | 113.0 ± 17.8a |
AS | 41.7 ± 3.10a | 22.4 ± 1.24b | 44.0 ± 4.18a | 21.8 ± 3.41b |
On the other hand, the available sulfur (AS) concentration was significantly lower in stony than in non-stony soils (p < 0.05) (Table 1). The pH of RZS in stony soils was significantly higher (p < 0.05) than that in non-stony soils, but there was no significant difference between the pH of non-stony and stony soils in BS samples. In addition, there was no significant difference between any of the physicochemical parameters between BS and RZS of the same soil type (p > 0.05). The mean tuber size of T. hemsleyanum grown in stony soils for two years was markedly larger than that in non-stony soils (Fig. S1 in Additional file 1). For example, the average length, diameter and weight of tubers grown in stony soils were 1.7-, 1.9- and 5.6-fold higher than those in non-stony soils (Fig. S1 in Additional file 1).
Differences In Soil Microbiota Between Non-stony And Stony Soils
For all 38 soil samples, a total of 1.63 million reads was obtained after quality filtering, representing 24,700 bacterial operational taxonomic units (OTUs) at 99% sequence similarity. The bacterial alpha-diversity indices, such as Chao1, observed species richness, phylogenetic diversity (PD), and Shannon diversity index, in stony soils were generally higher than those in non-stony soils, with both the BS and RZS samples showing significant higher values (p < 0.05) in stony soils, compared with the corresponding samples in non-stony soils.
In the soil samples from tuber surface (TS), the observed species richness and PD of stony soils were also markedly higher than the corresponding values of non-stony soils (Fig. S2A in Additional file 1). Principal coordinate analysis (PCoA) showed that the samples from non-stony and stony soils were clearly separated along the PCoA1 axis (Fig. S2B in Additional file 1). This distinction was further confirmed via the analysis of similarities (ANOSIM) and permutational multivariate analysis of variance (PERMANOVA) analyses, indicating that the bacterial communities of stony soils were significantly different from those of non-stony soils (p = 0.001) (Fig. S2B in Additional file 1). Moreover, the bacterial communities of BS, RZS, rhizosphere soil (RS), and TS in stony soils also showed significant differences (p < 0.05) with the corresponding communities in non-stony soils (Table 2).
Table 2
ANOSIM similarity analysis between non-stony and stony soils. ANOSIM similarity analysis based on Bray-Curtis dissimilarity for the comparison of bacterial community composition in the BS, RZS, RS and TS between non-stony and stony soils (permutation = 999). Values in bold indicated significant influences (p < 0.05).
| St-BS | Non-RZS | St-RZS | Non-RS | St-RS | Non-TS | St-TS |
R | P | R | P | R | P | R | P | R | P | R | P | R | P |
Non-BS | 1.00 | 0.030 | -0.02 | 0.462 | 1.00 | 0.007 | 0.56 | 0.018 | 1.00 | 0.011 | 0.43 | 0.04 | 0.57 | 0.009 |
St-BS | | | 0.66 | 0.011 | 0.08 | 0.198 | 0.74 | 0.017 | 0.71 | 0.016 | 0.69 | 0.034 | 0.21 | 0.098 |
Non-RZS | | | | | 0.60 | 0.019 | 0.18 | 0.111 | 0.72 | 0.016 | 0.06 | 0.245 | 0.41 | 0.014 |
St-RZS | | | | | | | 0.75 | 0.009 | 0.57 | 0.009 | 0.77 | 0.011 | 0.13 | 0.165 |
Non-RS | | | | | | | | | 0.49 | 0.039 | -0.20 | 0.994 | 0.38 | 0.018 |
St-RS | | | | | | | | | | | 0.51 | 0.014 | 0.19 | 0.089 |
Non-TS | | | | | | | | | | | | | 0.33 | 0.034 |
To explore whether the differences in microbial structure reflected the changes in bacterial community composition, we further analyzed the differences in taxonomic identity and relative abundance of the bacterial taxa for each soil. The most abundant phyla were Proteobacteria, Acidobacteria and Actinobacteria in both non-stony and stony soils (Fig. S3A in Additional file 1). However, the phylum Acidobacteria showed a lower relative abundance in stony than in non-stony soils, especially in the BS samples, whereas the relative abundance of the phylum Actinobacteria was higher in stony than in non-stony soils, especially in the RZS and TS samples. In addition, the relative abundances of phylum WPS-2 in BS, RZS, RS, and TS samples of stony soils were significantly lower (p < 0.01), and the abundances of phylum Rokubacteria in these soil samples were significantly higher (p < 0.01 in BS; p < 0.05 in the others), compared with those in non-stony soils (Fig. S3B in Additional file 1).
At the family level, the Methyloligellaceae and families belonging to the Actinobacteria in BS, RZS, RS, and TS samples of stony soils were significantly more abundant than those of non-stony soils (Welch’s t-test, p < 0.05, Storey false discovery rate (FDR)-corrected) (Fig. 1). Furthermore, the relative abundances of the Caulobacteraceae and families from WPS-2 in BS, RZS, and TS (Fig. 1A, C and D), and certain families from the phylum Gammaproteobacteria in RZS, RS, and TS of stony soils (Fig. 1B-D) were markedly lower than the corresponding values of non-stony soils (Welch’s t-test, p < 0.05, Storey FDR-corrected). Moreover, in BS samples, the families from the phyla Acidobacteria and Chloroflexi of stony soils had a significantly lower abundance than those of non-stony soils (Welch’s t-test, p < 0.05, Storey FDR-corrected) (Fig. 1A).
At the OTU level, there was a total of 110 OTUs, the relative abundances of which were more than 0.5% in at least one group (Fig. S4A in Additional file 1). Among them, there were 35 OTUs in BS samples which were significantly different between non-stony and stony soils, and, correspondingly, there were 16, 16 and 12 OTUs in RZS, RS and TS samples, respectively, based on three screening criteria (Fig. 2A and Additional file 2). The four samples shared five OTUs, including four OTUs belonging to the Alphaproteobacteria (Rhizobiales_Xanthobacteraceae_OTU2030, Rhizobiales_ Methyloligellaceae_OTU6287, Elsterales_OTU5022, and Elsterales _OTU5266), and one OTU from Actinobacteria (OTU7092) (Fig. 2B).
Twenty-one OTUs with significant differences in at least two groups included eight upregulated OTUs and 13 downregulated OTUs found in stony soils, compared with the corresponding OTUs in non-stony soils (Fig. 2B). These discriminatory OTUs can be identified as key drivers of variability to distinguish the non-stony and stony soils along the first axis, according to the non-metric multidimensional scaling analyses (Fig. S4B in Additional file 1). The total contribution rate of these 21 OTUs to the differences in bacterial communities in all samples between non-stony and stony soils was 34.95%, with corresponding values of 37.78, 34.98, 35.68 and 33.66% in BS, RZS, RS and TS samples, respectively (Fig. 2B). In order to demonstrate whether the difference in bacterial community in TS samples was caused by the stony soils, we performed a SourceTrack analysis (Fig. S5 in Additional file 1). The data showed that about 70% of the microbiota from RZS was sourced from BS, and the microbiota of RZS further changed the bacterial community compositions of RS and TS (Fig. S5 in Additional file 1).
Co-occurrence networks of non-stony and stony soils were constructed, based on Spearman’s rank correlation (Fig. 3). The results indicated that the co-occurrence network of stony soils had higher numbers of edges, average connectivity, graph density, and modularity, but had shorter average path length, compared with the network of non-stony soils (Table S1 in Additional file 1). The network of stony soils had more nodes that belonged to Alphaproteobacteria (31.0% vs 30.2%), Actinobacteria (17.8% vs 14.6%), Verrucomicrobia (5.0% vs 2.7%) and Rokubacteria (1.8% vs 0.4%), and possessed fewer nodes assigned to Acidobacteria (19.5% vs 23.1%), Chloroflexi (7.92% vs 8.3%), and Gammaproteobacteria (7.1% vs 10.1%) than that of non-stony soils (Fig. 3A-B). Netshift analysis revealed that the number of significant associations in the communities of stony soils was more than that of non-stony soils (277 vs 224 associations, respectively), whereas only 21 out of 522 OTU associations existed in both communities (Fig. 3C). Moreover, a total of 54 potential driver OTUs were related to the changes in microbiome composition from non-stony to stony soils (Fig. 3C). Among these OTUs, 18 out of 21 significantly different OTUs mentioned above were also found to play important roles in changing the networks’ structure (Fig. 2B, 3C and Additional file 3).
Linking the bacterial community with soil physicochemical parameters and tuber size
To understand the factors driving different bacterial communities, the effects of physicochemical parameters on the bacterial communities of BS, RZS, RS, and TS samples were determined by the Mantel tests and distance-based linear modeling (DISTLM) analysis (Table S2-3 in Additional file 1). The Mantel tests showed that environmental variables had a significant influence on the composition of bacterial communities of BS, RZS, and RS (p < 0.05), and the parameters pH, TN, OM, and AS all explained large variations of bacterial communities in BS, RZS and RS (Table S2 in Additional file 1). DISTLM analysis indicated that individual AS had the strongest correlation with the variation in BS (pseudo-F = 5.62, p = 0.002), RZS (pseudo-F = 3.33, p = 0.003) and RS (pseudo-F = 2.70, p = 0.003), with corresponding explained percentages of 48.4, 29.4, and 25.3%, respectively, followed by pH, with explained percentages of 41.7, 23.0 and 19.5%, respectively; while individual TN had the most significant influence on the variation of TS (pseudo-F = 1.87, p = 0.008), with an explained percentage of 18.9% (Table S3 in Additional file 1). The sequential model indicated that soil factors accounted for 48.4, 44.3, 39.2 and 31.3% of the total variations in the BS, RZS, RS, and TS samples, respectively (Table S3 in Additional file 1).
Pearson’s correlation analysis of 21 significantly discriminatory OTUs and soil physicochemical parameters indicated that most physicochemical parameters, especially AS, revealed a significant correlation with these discriminatory OTUs, with AS having the most associations with discriminatory OTUs, followed by pH and TN (Fig. 4). To evaluate whether these discriminatory OTUs were associated with tuber size of T. hemsleyanum, we further analyzed the correlations between these discriminatory OTUs and four tuber size parameters (Fig. 4). The results indicated that the OTUs belonging to Methylologenllaceae (OTU265 and OTU6287) and Acidothermaceae_Acidothermus (OTU1962) were significantly positively correlated with tuber size, whereas the OTUs associated with Xanthobacteraceae (OTU2030), Burkholderiaceae _Burkholderia (OTU2662) and Chloroflexi (OTU7503) had significantly negative correlations with tuber size (Fig. 4).
The transcriptome of tubers grown in non-stony and stony soils
To identify whether different soil types can affect the transcriptome of tubers of T. hemsleyanum, total RNAs from tubers produced in non-stony and stony soils were sequenced separately on an Illumina HiSeq™ platform. The results showed that a total of 44,817,581 clean reads were obtained, with a 50.4% GC content and a 98.7% Q20 score. The clean data were then assembled into 115,777 unigenes with an average size of 550 bp (Table S4 in Additional file 1). The assembled unigenes were annotated using BLAST, based on sequence similarity searches against nine different public databases, and a total of 70,194 (60.6%) unigenes were annotated to at least one significant match in all public databases (Table S5 in Additional file 1).
Correlation analysis showed that the gene expression patterns in tubers of T. hemsleyanum grown in the same soils were similar, but there was a marked difference in profile between non-stony and stony soils (Fig. S6 in Additional file 1). DESeq2 analysis indicated that a total of 3,853 differentially expressed genes (DEGs), consisting of 1,145 upregulated genes and 2,708 downregulated genes, were detected in tubers of T. hemsleyanum grown in stony soils, compared with those in tubers grown in non-stony soils (Fig. 5A and Additional file 4). To understand the biological function of these DEGs, Gene Ontology (GO) enrichment analysis was performed. The results revealed that a large number of significant GO terms were identified between non-stony and stony soils, and the top 30 upregulated terms and top 30 downregulated terms were selected for plotting (Fig. 5).
The biological processes of the top 30 upregulated terms were prominently involved in chloroplast fission, regulation of cellular respiration, xanthine catabolism, starch biosynthesis, abscisic acid (ABA) biosynthesis, and peptide biosynthesis (Fig. 5B). The top 30 downregulated terms were closely connected with GA- and ethylene-mediated signaling pathways, response to chitin and biotic stimuli, and shoot system development (Fig. 5C). Due to the large number of DEGs between non-stony and stony soils, we used Gene Set Variation Analysis (GSVA) to calculate the abundance of GO pathways related to each soil type. A total of 33 significantly different (p < 0.05) GO terms, mainly associated with auxin homeostasis, photosynthesis, secondary metabolism, response to various stresses, and organic substances, were found (Fig. S7 in Additional file 1).
The key taxa in TS potentially affected tuber size through association with various GO pathways
To determine the correlation between bacterial taxa and host pathways, a network based on the Spearman’s correlation analysis was constructed, using the Cytoscape software (Fig. 6). The results indicated that the key taxa in TS had a large number of associations with 25 differentially expressed GO terms. The significantly upregulated OTUs, belonging to the Methylologenllaceae (OTU265 and OTU6287), Beijerinckiaceae_Rhodobaastus (OTU7514) and Acidothermaceae_Acidothermus (OTU1962), showed very strong positive correlations with the pathways related to IAA and ABA biosynthesis, and photosynthesis, whereas the most downregulated OTUs were negatively associated with these pathways. On the other hand, the downregulated OTUs belonging to Burkholderiaceae_Burkholderia (OTU2662 and OTU7204), Acidobacteria (OTU2017), Elsterales (OTU5022 and OTU5014), Xanthobacteraceae (OTU2030), and Chloroflexi (OTU7503), had a number of positive correlations with the terms involved in response to GA, ethylene, various biotic stresses, and organic substances (Fig. 6).