Phenotypic differences between salt-tolerant and salt-sensitive genotypes under salinity stress
Compared with control, seed germination rates of QQ056 (salt-tolerant, ST) and 37TES (salt-sensitive, SS) under 300 mM NaCl treatment were both significantly decreased, but the descender of SS was obviously higher than ST (Fig. 1a and c). And the germination rate of ST was significantly higher than SS under salt stress. Likewise, the ability to resist salinity stress of SS seedlings was also weaker than ST. The two contrasting genotypes both showed reduced growth, leaf chlorosis, delicate seedlings and decreased leaf number in response to salt stress, but ST was less damaged than SS (Fig. 1b). Although no obvious differences in plant survival rate, plant height, root length, root fresh weight and root dry weight were observed between the two genotypes under the control condition, the reduction in plant survival rate and root dry weight of ST genotype were less significant than SS and these two indexes of ST were both significantly higher than SS under salt treatment (Fig. 1c). Taken together, ST genotype possesses better salt tolerance than SS genotype. Compared with control, soluble sugar content and proline content of the two genotypes were both significantly increased, the superoxide dismutase (SOD) activity was significantly rised in ST but slightly declined in SS (Fig. 1d). Except SOD activity, there were no significant differences in soluble sugar content and proline content of both genotypes under control condition. But all of the three physiological indexes of ST were significatly higher than SS under salt treatment.
RNA sequencing and transcriptome assembly
The RNA samples from the roots of 28-day-old ST and SS seedlings under the control and salt-treated (300 mM NaCl, 0.5, 2, and 24 h) conditions were used for Illumina Genome Analyzer deep sequencing. On average, 41,470,921 raw reads were generated, and 38,510,203 clean reads were obtained after cleaning and quality checks from each sample (Table S1). The average clean rate can be high as 92.86%, demonstrating the high quality of the sequencing results. Approximately 85.30% of the reads were mapped to the quinoa reference genome, 33,182,448 reads were mapped to unique regions and 2,192,951 reads were mapped to multiple regions (Table S1).
Identification and analysis of DEGs under salt treatment
Compared with control, 780/1212, 8885/5940, and 6975/6328 genes were differentially expressed in ST/SS at 0.5, 2, and 24 h, respectively (Fig. 2a). The total number of DEGs was greater in ST than in SS under salt stress. However, the number of up-regulated and down-regulated DEGs in ST and SS at different salt stress time points also existed differences. Among these DEGs, 521/633, 4549/3299, and 2913/2779 were up-regulated in ST/SS at 0.5, 2, and 24 h, respectively, while 259/579, 4336/2641, and 4062/3549 were correspondingly down-regulated (Fig. 2a). Venn-diagram analysis indicated that the salt-responsive genes were genotype specific and time specific, which may account for the differences in salt tolerance of ST and SS genotypes (Fig. 2b and c).
The KEGG pathway enrichment analysis among all the DEGs were performed to reveal several important salt-related pathways. These pathways mainly involved biosynthesis of secondary metabolites, amino acid metabolism, vitamin metabolism, lipid metabolism, plant hormone signal transduction, carotenoid biosynthesis, and carbohydrate metabolism (Table S2). Among these, biosynthesis of secondary metabolites, metabolic pathways, plant hormone signal transduction and alpha-Linolenic acid metabolism were significantly enriched in both genotypes at each salt stress time point (Fig. 3). The significant enrichment of plant hormone signal transduction revealed the importance of plant hormones in salt stress. Based on 271 background genes in plant hormone signal transduction in quinoa, 25/30, 133/106, and 108/82 DEGs were enriched in this pathway in ST/SS at 0.5, 2, and 24 h, respectively (Table S2). In addition, some pathways were preferably enriched in ST or SS at all salt stress time points, such as plant-pathogen interaction in ST, amino sugar and nucleotide sugar metabolism and phenylpropanoid biosynthesis in SS, suggesting the differences in response mechanisms to salt stress of the two genotypes (Fig. 3).
GO enrichment analysis was conducted to characterize the biological functions of the DEGs. The results indicated that several GO terms were significantly enriched in both genotypes at all salt stress time points, including response to oxidative stress, oxidation-reduction process, oxidoreductase activity, and antioxidant activity, which have been proved to be related to salt tolerance in plants (Table S3 and S4). At 0.5 h of 300mM NaCl stress, transcription factor activity, sequence specific DNA binding, regulation of RNA biosynthetic process, and regulation of gene expression were significantly enriched in both genotypes, while the GO terms response to oxidative stress, peroxidase avtivity, and antioxidant activity were significantly enriched at 2 h and 24 h of salt stress (Fig. 4). These suggested that the initial stage of response to salt stress mainly involved the transcription activation or expression of some important genes with molecular functions, with the extension of salt treatment time some biological processes were induced. Meanwhile, the important role of transcription factors in response to salt stress was also revealed.
Analysis of the DEGs exclusively in ST genotype in response to salt stress
A total of 308 and 404 DEGs were commonly identified at all sampling time points in ST and SS, respectively, whereas 191 and 287 DEGs were exclusively found in ST and SS under NaCl treatment, respectively (Fig. 2b). Further analysis showed that 87 genes were up-regulated in ST among the 191 DEGs above, including calcium-dependent protein kinase, ethylene-responsive transcription factors, heat shock factor protein, UDP-glycosyltransferase, transcription factors WRKY and MYB, which might be candidate genes for improved salt tolerance in ST (Table S5). The most representative GO categories amongst the 191 DEGs exclusively in ST were biological processes related to response to oxidative stress, response to stimulus and response to stress, molecular functions related to peroxidase activity, antioxidant activity and sequence-specific DNA binding, cellular components related to extracellular region and cell wall (Table S6). The significantly enriched KEGG pathways included phenylpropanoid biosynthesis, biosynthesis of secondary metabolites and metabolic pathways (Table S7). 11 genes were enriched in phenylpropanoid biosynthesis, most of which were related to peroxidases.
Common DEGs between ST and SS under salt stress
In spite of the exclusively expressed differential genes in each genotype, 117 DEGs were common to various stages of both genotypes (Fig. 2b). These genes (83 up- and 34 down-regulated) were constitutively active and were called core salt-responsive genes despite different salt tolerance levels in quinoa, including many peroxidases (PODs), protein phosphatase 2Cs (PP2Cs), zinc finger proteins (ZATs), F-box proteins and plasma membrane ATPase (Table S8). Especially, several transcription factors were also included, such as WRKY, NAC and MYB, suggesting their important roles in transcription regulation contributing to salt tolerance of quinoa.
In addition, the expression patterns of the constitutively active genes were analyzed in both genotypes under salt stress. As shown in Fig. 5, the temporal and spatial expression patterns of 117 core DEGs were presented by a heatmap. Although there was a big difference in salt tolerance between the two genotypes, these DEGs exhibited similar expression patterns, indicating that these genes were commonly induced in this species in response to salinity stress.
GO enrichment analysis of common DEGs of ST and SS in response to salt stress
Based on GO enrichment analysis, the 117 core genes were categorized into three GO ontologies and 32 terms (Fig. 6). Many genes were significantly enriched in biological processes concerned with polysaccharide metabolic process, cellular glucan metabolic process, cellular polysaccharide metabolic process, glucan metabolic process, cellular carbohydrate metabolic process, carbohydrate metabolic process, and oxidation-reduction process (Fig. 6 and Table S9). Furthermore, some GO categories were significantly enriched in molecular functions amongst the DEGs, including xyloglucan:xyloglucosyl transferase activity, hydrolase activity, hydrolyzing O-glycosyl compounds, hydrolase activity, acting on glycosyl bonds, transferase activity, transferring hexosyl groups, transferase activity, transferring glycosyl groups, and oxidoreductase activity (Fig. 6 and Table S9). Hence, the salt resistance process is very complicated in quinoa.
KEGG pathway enrichment analysis of common DEGs of ST and SS in response to salt stress
The core genes were aligned with the KEGG database and were assigned to 16 KEGG pathways (Table S10). The significantly enriched pathways were alpha-Linolenic acid metabolism, biosynthesis of secondary metabolites, phenylpropanoid biosynthesis, plant hormone signal transduction, metabolic pathways, and SNARE interactions in vesicular transport (Table S10). The metabolic pathways exhibited the largest number of DEGs, indicating many biochemical reactions were activated to enhance the adaptation to salt environment. Four DEGs were enriched in phenylpropanoid biosynthesis, three of which were identified as PODs.
Five differentially expressed transcripts were predicted to participate in plant hormone signal transduction, including two PP2Cs (AUR62012136, AUR62040699), one TIFY10A (AUR62009919), and two ABA receptor PYLs (AUR62012310, AUR62022950) (Table S8 and S10). PP2C and PYL are important participants in ABA signaling pathway, their involvements eventually cause stomatal closure to respond to salt stress (Fig. 7). The JA-responsive protein TIFY10A can immediately induce ubiquitin mediated proteolysis or indirectly promote monoterpenoid biosynthesis and indole alkaloid biosynthesis to respond to stimulus (Fig. 7). The PP2Cs and TIFY10A were positively regulated in both genotypes under 300 mM NaCl treatment, while the PYLs were negatively regulated, indicating their different roles in the adaptation to high salinity. Meanwhile, it also revealed the significance of plant hormones in resistance to salt stress in quinoa.
Validation of RNA-Seq data by qRT-PCR analysis
To further validate the expression data of RNA-Seq, qRT-PCR was performed on the same RNA samples originally used for next-generation sequencing. 21 DEGs common in both genotypes under salt stress were selected for qRT-PCR analysis, including six transcription factors (four NACs, one MYB and one WRKY), five plant hormone-related genes (two PP2Cs, one TIFY10A and one PYL) and other ten key genes (Fig. 8a). Among these, 13 were up-regulated and 8 were down-regulated (Table S8). The specific primers are listed in Table S11. To compare the expression data between RNA-Seq and qRT-PCR, the relative expression level was transformed to log2 Fold Change. The qRT-PCR results were strongly correlated with RNA-Seq data both in ST (R2 = 0.8005) (Fig. 8b) and SS (R2 = 0.7973) (Fig. 8c), demonstrating the reliability of the RNA-Seq expression profile in this study.