Identification of DEGs involved in PSNPs treatments
Differentially expressed genes (DEGs) in wheat shoots and roots were identified by comparing the expression level of each gene in PSNPs treatments and the control group (Fig. 1a; Fig. S1). As shown in Fig. 1c, the number of differentially up-regulated or down-regulated genes in the shoots and roots of wheat steadily increased with the increase of PSNPs exposure concentration. At the same time, the number of genes down-regulated was almost twice that of up-regulated genes in wheat tissues. The most interesting aspect of Fig. 1c was that regardless of the PSNPs exposure concentration, the number of down-regulated genes in the roots is always higher than that in the shoots. Additionally, the Venn diagram further illustrates the possible DEGs relations among different PSNPs treatments, where 1910 (7.83%) and 1167 (6.11%) DEGs were common to all four PSNPs concentrations in root and shoot, respectively (Fig. 1d). Taken together, these above results suggested that PSNPs did induce genetic changes in wheat as environmental xenobiotics and meanwhile down-expressed genes were the dominant part of differential genes involved in PSNPs treatments, especially in wheat root.
Characterization of co-expression networks related to PSNPs-induced traits
In order to figure out the candidate regulatory genes related to our desired traits (agronomic characters, mineral nutrients, and metabolites), all filtered transcripts from different wheat tissues were pooled together for subsequent WGCNA analysis. Accordingly, a total of eight modules were obtained in our study that were named in different colors (Fig. 2a,b). Among all the modules, the turquoise module (MEturquoise) contained the largest number of genes (5678), followed by MEblue (5286), MEbrown (4559), MEyellow (2703), MEgreen (1009), and MEred (349) (Fig. S3). It could be seen that the expression level of the module in wheat was highly tissue-specific, which was also reflected in the expression patterns of DEGs (Fig. 2b; Fig. S1). By way of example, the eigengenes in blue and green modules were over-expressed in the roots, while low-expressed in shoots. Similar phenomena were also found in MEturquoise and MEred. Moreover, the hierarchical cluster analysis result demonstrated that MEblue and MEgreen were provided with high homogeneity, and they appeared to be closely related (Fig. S2c).
Through correlation analysis, we identified four modules (MEblue, MEgreen, MEturquoise, MEred) that were significantly related to the target traits (SPAD, biomass, Pn, Gs, Tr, and Cu, Mn, Zn, Fe, Mg contents) (Fig. 2c). A remarkably strong correlation was also found between these four modules and metabolites (Fig. S5). In particular, agronomic characters (SPAD, biomass, Pn, Gs, Tr) and almost all metabolites shared a common correlation trend, positively related to MEblue and MEgreen while negatively correlated with MEturquoise and MEred. However, for mineral nutrients, the results changed and depended largely on the element types. Obviously, it was worth noting that the relevance between macronutrient (Mg) and micronutrients (Cu, Mn, Zn, Fe) was element-dependent (Fig. 2c). In terms of eigengene expression within modules, transcripts comprising Meblue and MEgreen displayed high expression in root but low expression in leaf. Conversely, the transcripts expression trends in MEturquoise and MEred were totally opposite (Fig. S4). Overall, these results indicated that wheat gene expression profile exposed to PSNPs was tissue-specific, and functioning genes that control target traits were basically clustered in four modules, including MEblue, MEgreen, MEturquoise, and MEred.
Hub genes and GO terms associated with PSNPs exposure
To investigate the possible functions of tracked modules, we divided these four modules into two large sets where MEblue and MEgreen were regarded as the first groups, and MEturquoise and MEred were the second group. In this context, the top 20 genes were screened according to the degree of connectivity within each module for better display of the selected module co-expression network. Results of gene interaction analysis showed that TraesCS3D02G236700 and TraesCS3D02G353100 were the key genes of the blue and green modules respectively (Fig. 3a). Of note, no specific network functions of the module were identified via the hubs approach. Therefore, the genes in MEblue and MEgreen were merged for further Gene Ontology (GO) enrichment analysis, and the top 20 GO terms were presented as well (Fig. 3b). Most of the significant GO terms enriched by these modules were mainly nucleoside metabolic and protein processing processes, including GTPase activity (75), protein folding (53), cofactor metabolic process (63), nucleoside phosphate metabolic process (56), and purine nucleoside triphosphate metabolic process (47), etc. Other than these, some GO terms involved in metal transport were also found such as ion transmembrane transport (45), metal ion transport (47), cation transmembrane transport (30), inorganic ion transmembrane transport (30), copper ion transport (5), ion antiporter activity (5) and metal ion transmembrane transporter activity (24) (Supporting material excel_1). More importantly, the same GO terms for metal ion transport were also enriched by DEGs (Fig. S5, 6). In the case of metal ion transport, the over-expression of these eigengenes primarily occurred in wheat roots when exposed to high concentrations of PSNPs, while only a small number of genes showed the differential expression after PSNPs treatments (Fig. 3c).
As for the negative modules, key genes in turquoise (TraesCS7A02G459700) and red module (TraesCS3A02G481500) were also screened out, which belong to the Glycoside hydrolase family 17 (Fig. 4a, b). Correspondingly, the results showed that cofactor metabolic and biosynthetic processes (75, 45), sulfur compound metabolic process (38), amino acid metabolic and biosynthetic processes (47,74), and metal cluster binding (44) were the significant GO function items (Fig. 4c). With the exception of these, it was worth noting that many other GO terms are related to the metabolism and biosynthesis of amino acids, including serine family amino acid, glutamine family amino acid, aspartate family amino acid, sulfur amino acid, cysteine, lysine, histidine, and amino acid activation (Supporting material excel_2). In more detail, we investigated the expression profile of these eigengenes in different wheat tissues. Interestingly, the expression of these genes in roots did not change much while there were significant expression alterations in leaves especially at 0.1 mg/L PSNPs (Fig. 4d). As a result, it was reasonable to speculate that these two modules are mainly responsible for the biosynthesis and metabolism of amino acids.
Enriched KEGG pathway in wheat exposed to PSNPs
In addition to GO enrichment analysis, the mainly enriched biological pathways in
wheat induced by PSNPs exposure were identified by two independent KEGG analyses (Table 1). In WGCNA, there were 49(10.3%) and 47 (8.8%) genes in Meblue-green and Meturquoise-red module respectively that were assigned to carbon metabolism, while the most significantly enriched pathway in Meturquoise-red was ribosome. Notably, a large portion of genes in MEblue-green was matched to other biological pathways such as biosynthesis of amino acids, oxidative phosphorylation, glycine, serine and threonine metabolism, and pentose phosphate pathway. For comparison, it was expected to distinguish differences in enriched pathways in wheat shoot and root using MetaDE analysis. However, most of the DEGs were mapped to similar KEGG terms with slightly different p-value, which included plant-pathogen interaction, plant hormone signal transduction, MAPK signaling pathway – plant and glycerophospholipid metabolism.
Validation of the RNA-sequencing results
In the present study, the expression of 10 genes was verified by an RT-qPCR approach to assess the performance characteristics of RNA sequencing results. As shown in Fig. 5, there was a significant correlation (p < 0.01) in the relative gene expression between RT-qPCR and RNA-Seq method, indicating that the transcriptome sequencing results are reliable.