Skip to main content
Download PDF

Time-series transcriptome comparison reveals the gene regulation network under salt stress in soybean (Glycine max) roots

BMC Plant Biology volume 22, Article number: 157 (2022) Cite this article

Abstract

Background

Soil salinity is a primary factor limiting soybean (Glycine max) productivity. Breeding soybean for tolerance to high salt conditions is therefore critical for increasing yield. To explore the molecular mechanism of soybean responses to salt stress, we performed a comparative transcriptome time-series analysis of root samples collected from two soybean cultivars with contrasting salt sensitivity.

Results

The salt-tolerant cultivar ‘Qi Huang No.34’ (QH34) showed more differential expression of genes than the salt-sensitive cultivar ‘Dong Nong No.50’ (DN50). We identified 17,477 genes responsive to salt stress, of which 6644 exhibited distinct expression differences between the two soybean cultivars. We constructed the corresponding co-expression network and performed Gene Ontology term and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analysis. The results suggested that phytohormone signaling, oxidoreduction, phenylpropanoid biosynthesis, the mitogen-activated protein kinase pathway and ribosome metabolism may play crucial roles in response to salt stress.

Conclusions

Our comparative analysis offers a comprehensive understanding of the genes involved in responding to salt stress and maintaining cell homeostasis in soybean. The regulatory gene networks constructed here also provide valuable molecular resources for future functional studies and breeding of soybean with improved tolerance to salinity.

Peer Review reports

Background

Crops experience various abiotic stresses during their growth. Sodium (Na+) promotes plant growth at low concentrations but can be detrimental in high salinity conditions. Of the current 230 million hectares (ha) of irrigated land, 45 million ha (19.5%) are estimated to be affected by high salt; similarly, 32 million ha (2.1%) of the 1500 million ha used for dryland agriculture are affected by varying degrees of salt stress [1]. Affected areas are expected to increase by 10% annually due to changing global climatic conditions, low precipitation, high surface evaporation, weathering of native rocks, land clearing, irrigation with saline water, and poor cultural practices, reaching 50% of the world’s arable land by 2050 [2, 3]. When subjected to salt stress, plants experience osmotic stress, ionic toxicity and complex secondary effects [4]. High Na+ concentration in the soil leads to high osmotic pressure of the soil, disrupts cellular ion homeostasis, and prevents water and nutrient uptake from the soil, thus affecting plant growth and reducing yield [5, 6]. In addition, salt stress is accompanied by the accumulation of reactive oxygen species (ROS), which act as a secondary stress by causing the peroxidation of membrane lipids and the destruction of cellular membrane structures and proteins. Plants have accordingly developed multiple tolerance mechanisms to cope with salt stress, including adjustment and maintenance of ion homeostasis in response to osmotic stress, restoration of osmotic balance, and other metabolic and structural adaptations [6, 7].

Soybean (Glycine max (L.) Merr.), which originated in China and was domesticated ~ 6000–9000 years ago, is a vital source of both protein and cooking oil, providing 59% of the world’s oilseed production and 69% of the daily vegetable protein consumed [8]. Salt stress affects almost all aspects of soybean growth and development: germination, vegetative and reproductive growth, nodulation, leaf size, plant height, root length, shoot and root dry weight, seed size and weight [3, 7, 9]. Genes associated with tolerance to salt stress might be used in breeding new soybean varieties with high salt tolerance. Several genes participate in salt responses in soybean. These genes include ion regulator genes such as HIGH-AFFINITY K+ TRANSPORTER 1;1 (GmHKT1;1), GmHKT1;4, Arabidopsis K+ TRANSPORTER 1 (GmAKT1), CATION DIFFUSION FACILITATOR 1 (GmCDF1), qNaCl3 (GmNcl), Na+/H+ ANTIPORTER 1 (GmNHX1), GmNHX5, SALT OVERLY SENSITIVE 1 (GmSOS1) and CHLORIDE CHANNEL 1 (GmCLC1), which play an important role in maintaining ion homeostasis under salt stress by regulating the transport and accumulation of Na+, potassium (K+), chloride (Cl−), and other ions [10,11,12,13,14,15,16]. Protein kinase genes such as Ser/Thr PROTEIN KINASE 4 (GmPKS4), NIMA-RELATED KINASE 1 (GmNEK1) and CALCIUM-DEPENDENT PROTEIN KINASE 3 (GmCDPK3) [17,18,19], as well as transcription factor genes such as GmNAC, GmWRKY, basic leucine-zipper (GmbZIP), HEAT-SHOCK FACTOR (GmHSF), GmMYB, PLANT HOMEODOMAIN (GmPHD), GmDREB and NUCLEAR FACTOR Y SUBUNIT A (GmNFYA) [20,21,22,23,24,25,26] also respond to salt stress. Besides, phytohormones, including auxins, gibberellins, cytokinins, abscisic acid, ethylene, salicylic acid, jasmonates, brassinosteroids and strigolactones, also play major roles in mediating plant salt stress [6, 27]. However, few genes have been demonstrated to affect soybean yields in saline field conditions; little is known about salt signaling as soybean responds and adapts to high salinity. This lack of knowledge has greatly inhibited attempts to enhance salt tolerance in soybean.

Comparative differential gene expression analysis between contrasting genotypes is helpful to identify candidate genes and underlying molecular mechanisms. In this study, we selected the salt-tolerant soybean cultivar ‘Qi Huang No.34’ (QH34), a high-protein inbred line with high yields whose genome shows abundant genetic diversity inherited from six Chinese soybean accessions [28], and the salt-sensitive soybean cultivar ‘Dong Nong No.50’ (DN50), which was bred from the Canadian small soybean variety ‘Electron’ and is characterized by a short growth period, short plants, and single pods with many seeds. Here, we used transcriptome sequencing (RNA-seq) to identify salt stress-responsive genes in soybean. The identification of candidate genes and their associated mechanisms will provide a new basis for breeding salt-tolerant and high-yielding soybean varieties and contribute to increasing soybean productivity.

Results

Characteristics of QH34 and DN50 under salt stress

We evaluated the salt sensitivity of 22 soybean accessions (Table S1). Their phenotypes are shown in Fig. S1, in which Fig. S1a shows the phenotypes before salt treatment, Fig. S1b and S1c show the phenotypes grown in half-strength Hoagland’s solution without or with added NaCl, respectively. The growth of all 22 accessions was inhibited upon salt exposure (Fig. S1), but cultivar QH34 exhibited superior salt tolerance compared to the other accessions, while DN50 showed the lowest survival rate (Fig. 1a, b; Fig. S1). After 6 d of salt treatment, the leaves of QH34 remained green, whereas most leaves from DN50 turned yellow or died (Fig. 1a; Fig. S1). When plants were removed from salt stress and returned to half-strength Hoagland’s solution for 3 d, 84.5% of QH34 plants survived, whereas only 17.7% of DN50 plants survived (Fig. 1b). Root growth was also substantially slower under salt stress for both accessions, although the effect of salt stress was more severe for DN50 than for QH34, with the elongation of salt-treated roots repressed by 51.3% for DN50, but only 25.4% for QH34, relative to non-treated control roots (Fig. 1c). We therefore selected QH34 as a salt-tolerant accession and DN50 as a salt-hypersensitive accession for comparative root transcriptome analysis.

Fig. 1
figure 1

Salt tolerance in the two soybean cultivars DN50 and QH34 and RNA-seq analysis. a, Phenotypic differences between QH34 and DN50 in normal and salt-treated growth conditions. Bars = 5 cm. b, Survival rate of QH34 and DN50 plants after salt treatment. c, Percentage of root elongation repression by salt treatment. d, Hierarchical clustering of RNA-seq samples. e, PCA of RNA-seq samples

Full size image

RNA-seq data processing

We sequenced 24 root RNA libraries with an average number of 25.5 million 150-bp paired-end reads, ranging from 20.7 million to 30.6 million (Table 1). The mean rate of reads passing quality control was 99.3%, with an average GC content of 43.0% and duplicates of 43.3% (Fig. S2a; Table 1). We observed no significant difference in the sequencing output of the 24 samples, when estimated from Log2-normalized counts per million (Fig. S2b). We mapped the resulting clean reads to the soybean William 82 av2 reference genome using tophat2, reaching mapping rates of 94.1–96.0%, among which the percentage of unique mapped reads varied from 89.1 to 93.0% (Table 1), indicating that the RNA-seq libraries are of high quality.

Table 1 Evaluation of RNA-seq data from the two soybean cultivars DN50 and QH34
Full size table

Before proceeding with the identification of differentially expressed genes (DEGs), we performed hierarchical clustering (Fig. 1d) and principal component analysis (PCA; Fig. 1e) to estimate the similarity between samples. As shown in Fig. 1d, biological replicates globally clustered together. We identified two distinct clusters: cluster I with QH34 samples subjected to salt stress for 4 h and 8 h, as well as control samples for both QH34 and DN50, and cluster II with all other salt-treated samples. Root samples exposed to salt for 4 h and 8 h appeared to be more similar to the corresponding control samples for QH34, indicating that QH34 is less affected by salt treatment, in agreement with the strong salt tolerance observed in this variety. By contrast, samples collected from QH34 after only 2 h of salt treatment grouped with all salt-treated DN50 samples. PCA also indicated a good reproducibility across the biological samples, as they each defined a single group per sample type. In addition, the two major principal components can well separate two materials as well as the same material with difference treatment (Fig. 1e). In summary, both hierarchical clustering and PCA reflected the good quality and repeatability of our samples.

Identification of DEGs

We proceeded to identify DEGs in two steps. First, we determined the salt-responsive genes by comparing salt-treated samples with the controls for each cultivar at three time points (2 h, 4 h and 8 h). We detected fewer DEGs for QH34 than for DN50 after 2 h of salt treatment, while the opposite was true after 4 h and 8 h (Fig. 2a). We identified 13,890 unique salt-responsive genes in QH34, compared to 12,098 unique salt-responsive genes in DN50 across all three time points. We also obtained lower average fold-change values for the DEGs in QH34 than for those in DN50 (Fig. 2b) across all time points, indicating that when the salt stress is perceived, more genes were responsive in QH34, but the fluctuation of the expression of genes was controlled at a smaller scale compared to DN50. Second, we identified DEGs between the two cultivars for each time point (Fig. 2c–e). Genes that showed significant differences in expression between the two cultivars but were not salt-responsive were not analyzed further, yielding a final list of 1231 DEGs after 2 h of salt treatments, 5038 DEGs after 4 h, and 2895 DEGs after 8 h. The smaller number of DEGs identified after 2 h suggested that similar gene networks are triggered by salt stress in the two cultivars at this early time point. Of the 7875 salt-responsive genes in QH34 and 10,270 salt-responsive genes in DN50 after 2 h of exposure to salt, only 1231 genes were differentially expressed between the cultivars, with 520 (42.2%) were shared by QH34 and DN50, while the remaining 492 (40.0%) and 219 (17.8%, Fig. 2c) were DN50- and QH34-specific DEGs, respectively. By contrast, the relative proportion of shared DEGs identified after 4 h and 8 h of salt treatment was lower than after 2 h, as 1203 of 5038 DEGs (23.9%) and 687 of 2895 DEGs (23.7%) were shared for the 4 h and 8 h salt treatments, respectively, suggestive of cultivar-specific responses to salt stress occurred (Fig. 2d, e). We also compared DEGs across the three time points and identified 233 genes showing significant expression differences between QH34 and DN50 at all three time points, which defined universal salt-responsive genes (uDEGs). Other DEGs, either specific to a single time point or shared by two time points, were considered as specific salt-responsive genes (sDEGs) (Fig. 2f). We obtained similar numbers of upregulated and downregulated genes after salt stress (Fig. 2g). The 233 uDEGs showed distinct expression patterns between QH34 and DN50 (Fig. 2f), most exhibiting consistent behavior across all time points (Fig. 2h).

Fig. 2
figure 2

Time-dependent comparisons of DEGs between QH34 and DN50. a, Number of salt-responsive genes identified at each time point. b, Boxplot showing the average changes of DEGs at each time point, shown as absolute Log2(fold change). c–e, Venn diagrams showing the extent of overlap between DEGs from QH34 and DN50 after 2 h (c), 4 h (d), or 8 h (e) of salt treatment. f, Venn diagram showing shared and unique DEGs across the three time points. g, Identification of upregulated and downregulated genes h, Heatmap showing the expression patterns of DEGs across the three salt treatment durations

Full size image

Classification of the expression patterns of the DEGs

We classified the expression patterns of the DEGs into eight categories by calculating their transcript abundance relative either to the control or to the previous time point (Methods), using QH34 as a reference. Overall, we observed similar expression patterns for DEGs between QH34 and DN50 (Fig. 3a). However, within each of the eight defined categories, an average of 25.2% of DEGs displayed a different expression pattern from that of its cohort, ranging from 19.5% (pattern 5) to 29.6% (pattern 6; Fig. 3b). Expression patterns 2 and 6 were present in the most DEGs, each accounting for 16% (Fig. 3c). Expression patterns of 2 and 6 were opposite to each other as well as patterns 3 and 8: expression of genes belonging to pattern 2 was upregulated at 2 h and 4 h and downregulated at 8 h, while expression of genes belonging to pattern 6 was first downregulated at 2 h and 4 h and then upregulated at 8 h after treatment. Genes belongs to pattern 3 were upregulated at 2 h and 8 h and downregulated at 4 h after salt treatment, while expression of genes belonging to pattern 8 were downregulated at 2 h and 8 h and upregulated at 4 h in QH34 (Fig. 3a). By further dissecting the expression pattern of DEGs belong to those four major groups showing different expression patterns between QH34 and DN50, we found, the majority of the genes were either upregulated or downregulated across all three time points in both cultivars (Fig. 3a, b). Comparing to the continuously up- or down-regulating the expression of a gene, the expression patterns 2,3, 6 and 8 are likely to maintain cell homeostasis. DEGs with differential expression between QH34 and DN50 are good candidates for a gene regulatory network and a starting point for studying the molecular mechanism contributing to salt tolerance in QH34.

Fig. 3
figure 3

Analysis of DEGs expression patterns. a, Using the expression patterns of DEGs in QH34 as reference (cyan), the expression of the corresponding genes in DN50 were plotted (red). b, Percentage of genes showing different expression patterns between QH34 and DN50. c, Distribution of each expression pattern type

Full size image

Construction of gene network, gene ontology (GO) and Kyoto encyclopedia of genes and genomes (KEGG) enrichment analysis

To facilitate interpreting the biological functions of the DEGs identified above and understand the molecular mechanism controlling salt tolerance in soybean, we clustered DEGs based on the Pearson’s correlation coefficient values calculated from their expressions and built undirected networks based on their patterns of co-expression for uDEGs and sDEGs independently. Of 233 uDEGs, 203 (87.1%) successfully clustered into three groups, containing 95 (40.8%), 83 (35.6%) and 25 (10.7%) DEGs each. For cluster I, the most enriched GO term was associated with the negative regulation of abscisic acid (ABA) signaling (Fig. 4a). Examination of the expression pattern of the key genes within this pathway, such as Glyma.02G250200, Glyma.11G222600, Glyma.02G302400, Glyma.13G090600 and Glyma.18G035000, we found that compared to the controls, the expression levels of those genes were both upregulation in QH34 and DN50, particularly at 2 h, although the induction of gene expression was much more modest in QH34 than in DN50 (Fig. S3). Since these genes are negative regulators of ABA signaling, we hypothesize that ABA levels may be lower in DN50 upon salt stress compared to QH34. The most enriched GO term in Cluster II was associated with oxidoreductase activity, responses to oxidative stress (Fig. 4b). The expression of this set of genes was higher in QH34 than in DN50 for control samples (Fig. S4). However, the induction of gene expression by salt stress was more pronounced in QH34 compared to DN50, such as observed for Glyma.16G055900. For genes whose expression was repressed by exposure to salt, the repression was weaker in QH34 compared to DN50, as with Glyma.13G306900 and Glyma.12G195600. Conversely, for genes whose expression was lower in QH34 than DN50 in the control samples, the decline in gene expression was less marked in QH34 than in DN50, such as observed for Glyma.08G321100 (Fig. S4). Cluster III mainly comprised genes involved in the negative regulation of proteolysis (Fig. 4c; S5). Other GO terms enriched for uDEGs were involved in cellular ion homeostasis, signal transduction, protein modification and carbohydrate metabolism, etc. (Fig. S6). KEGG enrichment of uDEGs indicated that six pathways are significantly enriched: glycolysis, phenylpropanoid biosynthesis, glutathione metabolism, carbon fixation in photosynthetic organism, plant hormone signal transduction and mitogen-activated protein kinase (MAPK) signaling pathway (Fig. S7).

Fig. 4
figure 4

GO enrichment analysis for uDEGs. a, Negative regulation of ABA-activated signaling. b, Oxidoreductase activity signaling. c, Negative regulation of proteolysis signaling

Full size image

We also performed a functional characterization of sDEGs at 2 h, 4 h and 8 h after treatment, resulting in gene networks with three, four and four clusters, respectively (Fig. S8). As shown in Fig. 5, we identified pathways specific to each time point as well as pathways shared across all three time points. Phenylpropanoid biosynthesis, plant hormone signal transduction and MAPK signaling pathway were shared by all three time points, although the underlying genes varied across the time points. As previously mentioned, we identified only 1231 DEGs at 2 h, and the most significant pathway enriched at this time point was plant hormone signal transduction, followed by MAPK signaling pathway (Fig. 5a). Ribosome metabolism became the predominant pathway at both 4 h and 8 h (Fig. 5b, c), prompting us to examine the expression pattern of genes involved in the most enriched pathway at each time point between QH34 and DN50 (Fig. S9, 10, 11 and 12). Of the genes involved in plant hormone signaling, we detected 20 (43.5%) responsive to auxin stimulus that appear to encode four auxin/indole-3-acetic acid (IAA) family proteins based on sequence similarity (Fig. S9). Genes similar to Arabidopsis IAA4 (Glyma.19G161000, Glyma.02G000500) and IAA26 (Glyma.07G015200) were upregulated upon salt stress in DN50, while IAA27 (Glyma.08G203100) was downregulated in QH34 under the same conditions. By contrast, a gene with similarity to AUXIN-RESPONSE FACTOR 5 (ARF5; Glyma.14G217700) was significantly induced in QH34 upon salt stress, while ARF9 (Glyma.07G134800) was strongly downregulated in DN50. Taken together, these results indicate that auxin contributes positively to salt tolerance in soybean during the tested period. Genes specifically involved in gibberellic acid (GA; Glyma.14G086600, Glyma.11G155100), salicylic acid (SA; Glyma.15G232000), ethylene (ET; Glyma.08G105000), jasmonic acid (JA; Glyma.10G031900, Glyma.14G086600) and brassinosteroids (BRs; Glyma.13G224300) were also more highly expressed in QH34 than in DN50 (Fig. S9). The expression of 73.9% (34 of 46) of genes involved in MAPK signaling pathway was induced by salt stress, with higher levels in QH34 than in DN50 for 80.4% (37 of 46; Fig. S10). Examination of genes involved in phenylpropanoid biosynthesis revealed genes encoding caffeic acid O-methyltransferase (CMOT; Glyma.06G295700, Glyma.12G109800, Glyma.20G003500) as being highly induced at 4 h and 8 h in the salt-tolerant cultivar QH34 (Fig. S11). We also observed significant expression differences for 133 genes involved in ribosome metabolism at 4 h and 8 h, as they were expressed at lower levels in the salt-tolerant cultivar QH34 than in DN50 (Fig. S12).

Fig. 5
figure 5

KEGG pathway analysis for sDEGs. a, b and c show pathway analysis at 2 h, 4 h and 8 h, respectively. 1, extracellular; 2, intracellular; 3, plasma membrane; 4, nucleus; 5, undefined

Full size image

Functional validation of selected DEGs

We validated by Real Time Quantitative PCR (RT-qPCR) the expression patterns of selected seven DEGs (Fig. 6; Table 2). Five genes (Glyma.06G202300, Glyma.01G187700, Glyma.08G055900, Glyma.03G117000 and Glyma.18G281700) showed basically consistent results by RNA-seq and RT-qPCR. We selected two of them (Glyma.06G202300 and Glyma.01G187700) with high expression levels in QH34 for functional validation using the soybean hairy root system in the salt-sensitive cultivar DN50 (Fig. 7a, b). A phylogenetic analysis of Glyma.01G187700 and orthologous Arabidopsis proteins (Fig. 7a) revealed that Glyma.01G187700 encodes caffeoyl-CoA O-methyltransferase 1 (GmCCoAOMT1). Under salt stress, Glyma.01G187700 was more highly expressed in QH34 than in DH50 (Fig. 6a, h). Overexpression of Glyma.01G187700 conferred higher tolerance to salt than the empty vector control, as evidenced by the greater root elongation observed in the presence of 100 mM NaCl but not in the absence of salt stress (Fig. 7c–f). We also characterized above-ground phenotypes of the transgenic plants. When treated with 200 mM NaCl for 12 d, the above-ground parts of OE plants also exhibited higher vitality than empty vector controls, as most leaves remained green, in sharp contrast to the yellowing or dying leaves of the controls (Fig. 7g).

Fig. 6
figure 6

Validated DEGs by RT-qPCR. a–g, Relative expression levels of DEGs Glyma.01G187700 (a), Glyma.06G202300 (b), Glyma.08G055900 (c), Glyma.02G093200 (d), Glyma.18G281700 (e), Glyma.03G117000 (f) and Glyma.13G173400 (g), as determined by RT-qPCR. h–n, FPKM values for DEGs Glyma.01G187700 (h), Glyma.06G202300 (i), Glyma.08G055900 (j), Glyma.02G093200 (k), Glyma.18G281700 (l), Glyma.03G117000 (m) and Glyma.13G173400 (n)

Full size image
Table 2 OE and RT-qPCR primers used for functional validation of DEGs
Full size table
Fig. 7
figure 7

Phylogenetic analysis and functional validation of GmCCoAOMT1 and GmCYP75B1. a and b, Phylogenetic analysis of Glyma.01G187700 and Glyma.06G202300 from soybean and Arabidopsis orthologous proteins, respectively. c, Relative expression of Glyma.01G187700 and Glyma.06G202300 in transgenic soybean hairy roots, as determined by RT-qPCR. d and e, Hairy root elongation between empty vector (EV) and overexpression (OE) plants in normal and high salt conditions. f, Phenotypes of hairy roots between EV and OE in normal and high salt conditions. g, Phenotypes of above-ground parts of EV and OE plants before (upper) and after (lower panel) salt treatment. NS, not significant. **p < 0.01. Bar = 5 cm

Full size image

Glyma.06G202300 encoded a member of the cytochrome P450 superfamily (GmCYP75B1), likely involved in oxidation–reduction process (Fig. 7b). The expression of Glyma.06G202300 was induced by salt stress and continuously increased at 4 h and 8 h after treatment in both QH34 and DN50, though with a lower amplitude in DN50 (Fig. 6b, i; S4). Overexpression of Glyma.06G202300 in transgenic soybean hairy roots increased its transcript levels over four-fold relative to the empty vector control, as measured by RT-qPCR (Fig. 7c). Root elongation were slightly, but not significantly, longer than empty vector controls when Glyma.06G202300 OE hairy roots were grown in control conditions consisting of half-strength Hoagland’s solution (Fig. 7d). After salt stress exposure (100 mM NaCl) for 3 d, root elongation is significantly inhibited, and the OE of Glyma.06G202300 was able to reduce the inhibition of root elongation caused by the salt stress (Fig. 7e, f). Likewise, the above-ground part of Glyma.06G202300 OE plants exhibited higher vitality compared to the empty vector controls after being treated with 200 mM NaCl for 12 d, with most leaves remaining green (Fig. 7g).

Discussion

Plant responses to salt stress are complex, as evidenced by the > 17,000 genes whose expression changed under salt stress over the short time period examined in the present study. Crosstalk between pathways often further complicates the dissection of the underlying molecular mechanism. Here, we performed GO term and KEGG enrichment analysis on clustered gene co-expression networks and focused on identifying pathways with the greatest differences between the salt-tolerant cultivar QH34 and the salt-sensitive cultivar DN50. This analysis revealed several candidate pathways with a role in salt responses, as detailed below.

Oxidoreduction

Oxidoreduction reactions involve the transfer of electrons between donor and acceptor molecules and are essential for basic life functions including photosynthesis and respiration. Peroxidase can detoxify H2O2 and respond to environmental stresses such as wounding, pathogen attacks and oxidative stress [29]. Here, the main pathway identified for uDEGs Cluster II was oxidoreductase activity (Fig. 4b). The genes in this cluster were more highly expressed in the salt-tolerant cultivar QH34 than in the salt-sensitive cultivar DN50, indicating their positive role in the salt tolerance of soybean (Fig. S4).

MAPK

MAPKs are activated by various biotic and abiotic stresses. Salt stress induces two well-characterized MAPKs activating signaling molecules, phosphatidic acid and ROS via nicotinamide adenine dinucleotide phosphate (NADPH)-oxidase [30]. Glycine max MAP kinase 1 (GMK1) is regulated by phosphatidic acid and H2O2 during salt stress [31]. The three most differentially expressed genes within the MAPK signaling pathway were MAPKK2 (Mitogen-activated protein kinase kinase 2), BAK1 (BRI1-associated receptor kinase 1) and RboH (Respiratory burst oxidase homolog). All three genes were significantly upregulated in the salt-tolerant cultivar QH34 and slightly induced in the salt-sensitive cultivar DN50. In Arabidopsis, MAPKK2 contributes to cold and salt signaling by regulating MAPK6 and MAPK4 [32]. Overexpression of poplar (Populus trichocarpa) MAPKK enhances salt tolerance in tobacco (Nicotiana tabacum) [33]. BAK1 is the co-receptor of the membrane receptor BRASSINOSTEROID INSENSITIVE 1 (BRI1); under salt stress conditions, BAK1 works with BRI1 to transduce the BR signal to BR SIGNALING KINASE 1 (BSK1) and activate the phosphatase BRI1 SUPPRESSOR 1 (BSU1). BSU1 then inhibits the kinase BR-INSENSITIVE 2 (BIN2) and promotes the translocation of the transcription factors BRASSINAZOLE-RESISTANT 1 (BZR1)/BRI1 EMS SUPPRESSOR 1 (BES1) to the nucleus to induce the expression of BR-responsive genes, enhancing salt tolerance [34]. The RboH identified here is homologous to AtrbohF, which functions in maintaining cellular Na+/K+ homeostasis under salt stress, as a mutation in AtrbohF caused imbalance of Na+/K+ homeostasis, resulting in salt sensitivity [35]. In wheat (Triticum aestivum), high Na+ concentrations induce NADPH oxidase-dependent ROS generation, which elevates Ca2+ levels in roots [36]. The coordination of the MAPK pathway with plant hormones may play a key role in perceiving and transmitting external signals and maintain cellular homeostasis.

Phenylpropanoid biosynthesis genes

Phenylpropanoids are a group of plant secondary metabolites derived from phenylalanine that play crucial roles in the plant life cycle. Phenylalanine is first converted to cinnamic acid by deamination, and the resulting product is hydroxylated and methylated to generate coumaric acid and other acids with a phenylpropane unit. Aldehydes and alcohols can be produced when the reduction of CoA-activated carboxyl groups occurs. The alcohols are called monolignols and are the building blocks for lignin, a major component of the cell wall. Phenylpropanoid biosynthesis is important under salt stress [37]. In the present study, 76.9% of DEGs involved in phenylpropanoid biosynthesis were induced by salt stress, while 84.6% of the same DEGs showed on average a higher expression level in the salt-tolerant cultivar QH34 under salt stress compared to DN50, suggesting their positive roles in salt tolerance. The most upregulated genes (Glyma.06G295700, Glyma.12G109800, Glyma.20G003500) were annotated as encoding CMOT, orthologs to Arabidopsis ATOMT1 (At5g54160). The cell walls of omt1 mutants were more susceptible to enzymatic hydrolysis and displayed a higher digestibility compared to the wild type [38]. In the case of maize (Zea mays), mutation in the gene also led to the improved digestibility and decreased lignin levels of the cell wall [39]. Transgenic Arabidopsis plants overexpressing CMOT from Carex rigescens), a stress-tolerant grass, not only exhibited enhanced salt stress tolerance but also produced more lateral roots and accumulated more proline and chlorophyll contents than the wild type [40]. The high induction of CMOT expression observed in QH34 under salt stress may contribute to the outstanding salt tolerance of this cultivar.

Plant hormones and salt tolerance

Plant hormones regulate plant growth and development, as well as responses to abiotic and biotic stresses [6, 27]. In the present study, we determined that the ABA as the most enriched plant hormone signaling pathway for uDEGs while auxin as the most significant pathway enriched in sDEGs. The roles of ABA in salt stress have been repeatedly reported. ABA is the most efficient plant hormone for lowering the concentrations of Na+ and Cl− and the Na+/K+ ratio and for increasing K+ and Ca2+ concentrations, proline accumulation, and soluble sugar contents [41]. Hormonal priming with ABA (50 ppm) increases wheat capacity to cope with salinity [42]. Similarly, ABA-treated wheat plants exhibit higher salt tolerance and a significant decrease in Na+ concentrations in flag leaves [43]. A salt-tolerant variety of rice (Oryza sativa) displays a higher ability than a sensitive cultivar to produce ABA, underscoring the positive contribution of ABA to reducing stress effects [44]. In agreement, endogenous ABA levels increase in the rice root subjected to salt stress [45]. In soybean, ABA content in leaves similarly increases with exposure to elevated salt stress [46]. In the present study, we established that genes involved in negative regulation of ABA biosynthesis are upregulated in DN50, but only slightly upregulated or downregulated in QH34 (Fig. S3). We hypothesize that ABA content is higher in QH34 than in DN50, which may contribute to the greater salt tolerance of QH34.

Auxin, which plays major roles in plant growth, controls root growth and the proliferation of lateral roots. Local auxin minima serve as signals that trigger the transition from cell division to cell differentiation in the Arabidopsis root [47]. Upon salt stress, we observed that negative regulators of auxin signaling are either significant upregulated or downregulated to a lesser extent in the salt-sensitive cultivar DN50, while positive regulators of auxin signaling were mainly upregulated in the salt-tolerant cultivar QH34 (Fig. S9). The exogenous application of IAA enhances the ability of maize plants to respond to abiotic stress [48]. Although the role of auxin in salt stress is likely to be complex by controlling the balance between growth and stress response and crosstalk with other phytohormones, auxin positively contributed to the salt resistance observed in QH34 during the period examined here.

Other phytohormones such as GA, SA, JA, ET and BRs may also play important roles in adjusting plant growth to survive in high salt conditions, and complex crosstalk between phytohormones was anticipated. However, as genes related to these pathways were less enriched and showed no significant differences in expression between QH34 and DN50, they are not further discussed here.

The role of ribosome metabolism in salt tolerance

Ribosome biogenesis is a central process in any cell. In rapidly growing cells, high levels of ribosomal proteins (RPs) are translated. Maintaining cell functions requires a tight coordination between ribosomal RNAs and RPs, and the alteration of any step along the ribosome assembly process may negatively influence cell growth and cause proteotoxic stress [49]. Our comparison of the root transcriptome of the salt-tolerant cultivar QH34 and the salt-sensitive cultivar DN50 indicated that most genes participating in ribosome assembly and metabolism are strongly repressed at 4 h and 8 h after salt treatment in QH34, and the expression depression of some essential structural constituent of ribosome, such as Glyma.12G162300 (encoding 30S ribosomal protein S20), Glyma.11G201700 (encoding 50S ribosomal protein L17) could be observed at as early as 2 h in QH34 (Fig. S12). In comparison, the same genes remained unchanged or slightly downregulated at 2 h in the salt-sensitive cultivar DN50. Our results suggest that lower ribosome activity may be crucial during the early stages of cell response to salt stress. In the common ice plant (Mesembryanthemum crystallinum), which can tolerate Na+ concentrations exceeding that found in sea water and complete its life cycle, the ribosome-inactivating protein gene RIP1 is highly induced upon salt stress [50]. The capacity of the cell to mount a timely and appropriate response to stress is critical for determining cell fate. Typically, the stress responses first allow cells to adapt and recover from the specific stress; if the stimulatory insult cannot be resolved, cell death signaling pathways will be initiated [51]. The timely adjustment of ribosome activity observed in QH34 may help cells better adapt to salt stress or act as a signal for activating downstream pathways, as observed in a process referred to as the ribotoxic stress response [52].

Conclusions

Our root transcriptome comparison analysis between the salt-tolerant cultivar QH34 and the salt-sensitive cultivar DN50 at three time points over the course of salt stress showed that more genes display an active response to salt stress in QH34 but with smaller amplitudes relative to DN50. We identified 17,477 unique salt stress responsive genes across three time points, of which 6644 genes differentiate QH34 and DN50 at least at one of the time points tested, suggesting that salt tolerance in soybean is complex. Further clustering, gene network construction and GO/KEGG enrichment analysis of DEGs suggested that during the time period observed, genes involved in plant hormone signaling, oxidoreduction, phenylpropanoid biosynthesis, MAPK signaling and ribosome metabolism may play crucial roles in response to salt stress. We also validated one of the DEGs, as evidenced by the observed salt tolerance in DN50 resulting from its overexpression. In summary, the findings reported here advance our understanding of the molecular mechanisms regulating salt tolerance in soybean. The identified genes will be useful in breeding new soybean varieties with improved salt tolerance.

Methods

Plant materials and application of salt stress

The salt-tolerant cultivar Qi Huang No.34 (QH34) and the salt-sensitive cultivar Dong Nong No.50 (DN50) were selected for root deep transcriptome sequencing. Plants were grown in a growth chamber under 16 h of light and 8 h of dark (28 °C/20 °C, 50% relative humidity). Four days after sowing (vermiculite: nutritional soil = 1:1), seedlings were transferred to hydroponic growth conditions in half-strength Hoagland’s solution (pH 6.0), which was renewed every 3 d until reaching the VC growth stage (unifoliolate leaves are fully developed and the first trifoliolate leaf appears). Twenty plants selected from each cultivar at similar growth status were either subjected to 150 mM NaCl in half-strength Hoagland’s solution (pH 6.0) or maintained in half-strength Hoagland’s solution without NaCl as the control for 6 d. The roots were cleaned to remove traces of NaCl and the plants were then returned to half-strength Hoagland’s solution for 3 d, at which point the survival rate was calculated. Test of significance (p ≤ 0.05) was performed in R using Student’s t-test.

RNA extraction, library preparation and transcriptome sequencing

For each cultivar, root composite samples from ten plant were harvested 2 h, 4 h and 8 h after salt treatment or from the control plants (as described above), three biological replicates. Total RNA was extracted using the RNeasy Plant Mini Kit (QIAgene) according to the manufacturer’s instructions. The quality of RNA was evaluated on 1% (w/v) agarose gels; the integrity of the extracted RNA was then assessed using an Agilent Bioanalyzer 2100 system (Agilent Technologies) before library construction. High-quality RNA was used for RNA-seq library using the Truseq RNA Library Prep Kit (Illumina). RNA libraries were sequenced at Novogene (Beijing, China) to generate paired-end 150-bp reads on a HiSeq2000 platform.

Data processing and identification of DEGs

Raw data in fastq format were first quality control using Trimmomatic (ILLUMINACLIP: Adapter_seq.fa:2:30:10) [53] to remove sequencing adaptors, low-quality bases (Phred score ≤ 20), and reads containing ploy-N. Quality-trimmed data were further evaluated using fastQC to calculate Q20, Q30 and GC contents of the clean data. For the identification of DEGs, the soybean (Glycine max) William 82 av2 reference genome and its corresponding annotation file were retrieved from soybase (https://www.soybase.org/). Genome index was built using Bowtie2 [54]. Paired-end reads having passed quality controls were aligned to the soybean reference genome using Tophat2 [55] with default parameters. The resulting bam files were used for reads duplication rate estimation with R package DupRadar. Median-normalized and Log2-transformed data was used for hierarchical clustering and principal component analysis to estimate the similarity between samples. FeatureCounts was used to count the number of reads mapping to each gene. Identification of DEGs between treatments was performed with the R package DESeq2 [56]. Genes showing at least two-fold expression changes (Log2(fold change) ≥1 or ≤ − 1) and with an adjusted p value ≤0.01 were considered as differentially expressed. To classify these DEGs and construct the gene network contributing to salt tolerance in QH34, salt-responsive genes were identified by comparing salt-treated samples with their corresponding controls for each cultivar. Genes showing significant expression changes between the control and the salt-treated plants were considered as salt-responsive genes. Only salt-responsive DEGs were further analyzed between QH34 and DN50 to identify genes contributing to salt tolerance in QH34.

Analysis of DEGs expression patterns

The expression patterns of DEGs induced by salt stress in each cultivar were calculated separately. The expression changes for a gene were calculated according to the ratio of fragments per kb of exon model per million mapped fragments (FPKMs) for a gene after salt treatment (2 h, 4 h, or 8 h) relative to FPKMs for the matching control (CK) or relative to another salt treatment time that less it (2 h or 4 h). Depending on whether it is positive (upregulated) or negative (downregulated) relative to its expression at the preceding time point, it was classified into 8 types.

Gene network construction, GO functional enrichment and KEGG pathway analysis

To construct gene co-expression networks and study the coordination of gene expression under salt stress, a raw gene count table generated from FeatureCounts was normalized by median normalization using the EBSeq R package [57], then Log2 transformation was applied to the normalized data and genes with low count were removed from downstream analysis (threshold was set based on the normalized reads counts< 10). Calculation of correlation coefficients was performed with the R package psych (method = “ Pearson ”). The resulting correlation tables were further filtered based on absolute correlation (> 0.9) and adjusted p value (< 0.01). Network statistics were calculated using the R package igraph [58]. The visualization and clustering gene network were performed using Cytoscape. Briefly, the calculated correlation table and the network statistics tables were imported into Cytoscape, then the loaded data was further analyzed using clusterMaker [59] (cluster method = “community clustering (Glay)”) with default parameters. The clustered gene network was visualized using Degree Sorted Circle Layout to identify the hub genes (the top most connected genes in each cluster). GO functional enrichment and KEGG [60] pathway analysis were conducted in ClueGO; gene IDs were converted to Ensembl gene IDs (Glycine max [3847]). For GO functional enrichment, we focused on Biological Process. For both GO and KEGG enrichment analysis, the minimum number of genes was set to five for a GO term to be significantly enriched, and groups with > 50% overlap were merged. GO term/pathway network connectivity (kappa score) was set to Medium (0.5). DEGs were classified into two groups: one group containing DEGs annotated as salt-responsive shared by both QH34 and DN50, and the other containing time point-specific genes for QH34. For both groups, GO function and KEGG pathway enrichment analysis were independently performed for each gene network cluster.

RT-qPCR analysis

Total RNA was extracted using the Trizol reagent. The cDNA was synthesized following the instructions of the HiScript® II Q RT SuperMix for qPCR (+gDNA wiper) synthesis kit (Vazyme, Nanjing, China). qPCR was performed on a Step One Plus™ Real-Time PCR system (Applied Biosystems, USA) using ChamQ Universal SYBR qPCR Master Mix (Vazyme, Nanjing, China). Three independents biological replicates were analyzed. Samples used for RNA extraction were collected from plant roots treated as mentioned above.

Overexpression vector construction and plant transformation

The seven genes were filtered for RT-qPCR validation based on their functional annotations, expression patterns, respective pathways and the data from literature reports [34, 61,62,63,64]. Two genes of them were filtered for functional validation, and the primers (Table 2) used for target gene amplification were designed on the Vazyme website (https://crm.vazyme.com/cetool/singlefragment.html) and the coding sequences (Williams 82) were amplified. The resulting PCR products were ligated into the 35S-pBI121 (kanamycin resistance) vector and transformed into Escherichia coli DH5α by freeze-thaw method. Plasmids with the correct insert and free of mutations were introduced into Agrobacterium (Agrobacterium rhizogenes) strain K599 (streptomycin-resistant) using the freeze-thaw method. K599 colonies with the target construct were used to inoculate 10 mL YEP liquid medium and incubated at 28 °C for 10 h with shaking (220 rpm). After 10 h, 100 μL bacterial culture was spread onto a YEP solid plate and cultured in an incubator (28 °C) for 2 d for subsequent hairy root experiments.

Validation of gene function through soybean hairy root system

We generated transgenic hairy roots overexpressing the target gene (OE) and evaluated their salt sensitivity by growth in half-strength Hoagland’s solution with or without added NaCl, side by side with control hairy roots transformed with empty vector (EV). The salt-sensitive cultivar DN50 was used for functional validation by overexpressing the candidate genes. DN50 was sown in soil (a 1:1 mix of vermiculite and nutritional soil). Six days after sowing, the tip of a 1-mL medical syringe was used to puncture four wounds 0.5–1 cm below the cotyledons in a cross pattern. Then, cultures from the relevant Agrobacterium colonies were smeared onto the four wounds with pearl wool. The seedlings were covered with a clear plastic lid to maintain high humidity. Six days later, the plastic cover was removed and the infected area was covered with vermiculite to maintain humidity. Then, 14 d later, the taproot was cut off when the hairy root was ~ 5–10 cm in length and cultured in half-strength Hoagland’s solution (pH 6.0) for 3 d to adapt environment, and consistent seedlings were selected for salinity treatment of 3 d with 100 mM NaCl; controls received no NaCl. The root length before and after salt treatment was measured and the root relative elongation was calculated to evaluate salt tolerance. In addition, about 35 d after Agrobacterium infection, overexpressing hairy root plants and the control were transplanted into mixed soil. After 3 d, the hairy root plants were irrigated with half-strength Hoagland’s solution containing 200 mM NaCl, or half-strength Hoagland’s solution without NaCl as the control group. Each plant was irrigated with 50 mL solution every 2 d. After 12 d of salt treatment, the above-ground phenotypes were determined.

Availability of data and materials

The transcriptome data was submitted to NCBI (https://www.ncbi. nlm.nih.gov/geo/), and could be accessed using project ID/accession ID PRJNA766706. 22 soybean germplasm resources have been authorized.

Abbreviations

DN50:

Dong Nong No.50

QH34:

Qi Huang No.34

DEGs:

Differentially expressed genes

PCA:

Principal component analysis

uDEGs:

Universal salt-responsive genes

sDEGs:

Specific salt-responsive genes

MAPK:

Mitogen-activated protein kinase

ROS:

Reactive oxygen species

NADPH:

Nicotinamide adenine dinucleotide phosphate

RPs:

Ribosomal proteins

ABA:

Abscisic acid

MAPKK:

Mitogen-activated protein kinase kinase

GA:

Gibberellic acid

SA:

Salicylic acid

ET:

Ethylene

JA:

Jasmonic acid

BRs:

Brassinosteroids

BRI1:

BR-insensitive 1

BAK1:

BRI1-associated receptor kinase 1

RboH:

Respiratory burst oxidase homolog

BSK1:

BR signaling kinase 1

BSU1:

BRI1 suppressor 1

CMOT:

Caffeic acid O-methyltransferase

GO:

Gene ontology

KEGG:

Kyoto encyclopedia of genes and genomes

RNA-seq:

Transcriptome sequencing

RT-qPCR:

Real time quantitative PCR

OE:

Overexpression

EV:

Empty voter

FPKMs:

Fragments per kb of exon model per million mapped fragments

References

  1. Kumar N, Soren KR, Bharadwaj C, Sneha PPR, Shrivastava AK, Pal M, et al. Genome-wide transcriptome analysis and physiological variation modulates gene regulatory networks acclimating salinity tolerance in chickpea. Environ Exp Bot. 2021;187:104478.

    Article  CAS  Google Scholar 

  2. Asif MA, Schilling RK, Tilbrook J, Brien C, Dowling K, Rabie H, et al. Mapping of novel salt tolerance QTL in an Excalibur × kukri doubled haploid wheat population. Theor Appl Genet. 2018;131(10):2179–96.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Shrivastava P, Kumar R. Soil salinity: a serious environmental issue and plant growth promoting bacteria as one of the tools for its alleviation. Saudi J Biol Sci. 2015;22(2):123–31.

    Article  CAS  PubMed  Google Scholar 

  4. Zhu JK. Abiotic stress signaling and responses in plants. Cell. 2016;167(2):313–24.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Hasegawa PM, Bressan RA, Zhu JK, Bohnert HJ. Plant cellular and molecular responses to high salinity. Annu Rev Plant Physiol Plant Mol Biol. 2000;51:463–99.

    Article  CAS  PubMed  Google Scholar 

  6. Yu ZP, Duan XB, Luo L, Dai SJ, Ding ZJ, Xia GM. How plant hormones mediate salt stress responses. Trends Plant Sci. 2020;25(11):1117–30.

    Article  CAS  PubMed  Google Scholar 

  7. Phang TH, Shao GH, Lam HM. Salt tolerance in soybean. J Integr Plant Biol. 2008;50(10):1196–212.

    Article  CAS  PubMed  Google Scholar 

  8. Sedivy EJ, Wu FQ, Hanzawa Y. Soybean domestication: the origin, genetic architecture and molecular bases. New Phytol. 2017;214(2):539–53.

    Article  PubMed  Google Scholar 

  9. Papiernik SK, Grieve CM, Lesch SM, Yates SR. Effects of salinity, imazethapyr, and chlorimuron application on soybean growth and yield. Commun Soil Sci Plant Anal. 2005;36:951–67.

    Article  CAS  Google Scholar 

  10. Chen HT, He H, Yu DY. Overexpression of a novel soybean gene modulating Na+ and K+ transport enhances salt tolerance in transgenic tobacco plants. Physiol Plant. 2011;141(1):11–8.

    Article  CAS  PubMed  Google Scholar 

  11. Chen HT, Chen X, Gu HP, Wu BY, Zhang HM, Yuan XX, et al. GmHKT1;4, a novel soybean gene regulating Na+/K+ ratio in roots enhances salt tolerance in transgenic plants. Plant Growth Regul. 2014;73:299–308.

    Article  CAS  Google Scholar 

  12. Wang XS, Zhao JL, Fang QW, Chang XC, Sun MY, Li WB, et al. GmAKT1 is involved in K+ uptake and Na+/K+ homeostasis in Arabidopsis and soybean plants. Plant Sci. 2021;304:110736.

    Article  CAS  PubMed  Google Scholar 

  13. Zhang W, Liao XL, Cui YM, Ma WY, Zhang XN, Du HY, et al. A cation diffusion facilitator, GmCDF1, negatively regulates salt tolerance in soybean. Plos Genet. 2019;15(1):e1007798.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Do DT, Chen HT, Htt V, Hamwieh A, Yamada T, Tadashi S, et al. Ncl synchronously regulates Na+, K+ and cl− in soybean and greatly increases the grain yield in saline field conditions. Sci Rep. 2016;6:19147.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Sun TJ, Ma N, Wang CQ, Fan HF, Wang MX, Zhang J, et al. A golgi-localized sodium/hydrogen exchanger positively regulates salt tolerance by maintaining higher K+/Na+ ratio in soybean. Front Plant Sci. 2021;12:638340.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Li Francisca WY, Wong FL, Tsai SN, Phang TH, Shao GH, Lam HM. Tonoplast-located GmCLC1 and GmNHX1 from soybean enhance NaCl tolerance in transgenic bright yellow (BY)-2 cells. Plant Cell Environ. 2006;29(6):1122–37.

    Article  Google Scholar 

  17. Ketehouli T, Zhou YG, Dai SY, Carther KFI, Sun DQ, Li Y, et al. A soybean calcineurin B-like protein-interacting protein kinase, GmPKS4, regulates plant responses to salt and alkali stresses. J Plant Physiol. 2021;256:153331.

    Article  CAS  PubMed  Google Scholar 

  18. Pan WJ, Tao JJ, Cheng T, Shen M, Ma JB, Zhang WK, et al. Soybean NIMA-related Kinase1 promotes plant growth and improves salt and cold tolerance. Plant Cell Physiol. 2017;58(7):1268–78.

    Article  CAS  PubMed  Google Scholar 

  19. Wang D, Liu YX, Yu Q, Zhao SP, Zhao JY, Ru JN, et al. Functional analysis of the soybean GmCDPK3 gene responding to drought and salt stresses. Int J Mol Sci. 2019;20(23):5909.

    Article  CAS  PubMed Central  Google Scholar 

  20. Bian XH, Li W, Niu CF, Wei W, Hu Y, Han JQ, et al. A class B heat shock factor selected for during soybean domestication contributes to salt tolerance by promoting flavonoid biosynthesis. New Phytol. 2020;225(1):268–83.

    Article  CAS  PubMed  Google Scholar 

  21. Wang F, Chen HW, Li QT, Wei W, Li W, Zhang WK, et al. GmWRKY27 interacts with GmMYB174 to reduce expression of GmNAC29 for stress tolerance in soybean plants. Plant J. 2015;83(2):224–36.

    Article  CAS  PubMed  Google Scholar 

  22. Yang Y, Yu TF, Ma J, Chen J, Zhou YB, Chen M, et al. The soybean bZIP transcription factorgene GmbZIP2 confers drought and salt resistances in transgenic plants. Int J Mol Sci. 2020;21(2):670.

    Article  CAS  PubMed Central  Google Scholar 

  23. Du YT, Zhao MJ, Wang CT, Gao Y, Wang YX, Liu YW, et al. Identification and characterization of GmMYB118 responses to drought and salt stress. BMC Plant Biol. 2018;18:320.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Wei W, Huang J, Hao YJ, Zou HF, Wang HW, Zhao JY, et al. Soybean GmPHD-type transcription regulators improve stress tolerance in transgenic arabidopsis plants. Plos One. 2009;4(9):e7209.

    Article  PubMed  PubMed Central  Google Scholar 

  25. Nguyen QH, Vu LTK, Nguyen LTN, Pham NTT, Nguyen YTH, Le SV, et al. Overexpression of the GmDREB6 gene enhances proline accumulation and salt tolerance in genetically modified soybean plants. Sci Rep. 2019;9:19663.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Lu L, Wei W, Tao JJ, Lu X, Bian XH, Hu Y, et al. Nuclear factor Y subunit GmNFYA competes with GmHDA13 for interaction with GmFVE to positively regulate salt tolerance in soybean. Plant Biotechnol J. 2021;19(11):2362–79.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Verma V, Ravindran P, Kumar PP. Plant hormone-mediated regulation of stress responses. BMC Plant Biol. 2016;16:86.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Huang WX, Hou JJ, Hu Q, An J, Zhang YW, Han Q, et al. Pedigree-based genetic dissection of quantitative loci for seed quality and yield characters in improved soybean. Mol Breeding. 2021;41:14.

    Article  CAS  Google Scholar 

  29. Fones H, Preston GM. Reactive oxygen and oxidative stress tolerance in plant pathogenic Pseudomonas. FEMS Microbiol Lett. 2012;327(1):1–8.

    Article  CAS  PubMed  Google Scholar 

  30. Jalmi SK, Sinha AK. ROS mediated MAPK signaling in abiotic and biotic stress-striking similarities and differences. Front Plant Sci. 2015;6:769.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Im JH, Lee H, Kim J, Kim HB, An CS. Soybean MAPK, GMK1 is dually regulated by phosphatidic acid and hydrogen peroxide and translocated to nucleus during salt stress. Mol Cells. 2012;34(3):271–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Teige M, Scheikl E, Eulgem T, Doczi R, Ichimura K, Shinozaki K, et al. The MKK2 pathway mediates cold and salt stress signaling in Arabidopsis. Mol Cell. 2004;15(1):141–52.

    Article  CAS  PubMed  Google Scholar 

  33. Yang CJ, Wang RN, Gou LZ, Si YC, Guan QJ. Overexpression of Populus trichocarpa mitogen-activated protein kinase kinase4 enhances salt tolerance in tobacco. Int J Mol Sci. 2017;18(10):2090.

    Article  PubMed Central  Google Scholar 

  34. Zhao SS, Zhang QK, Liu MY, Zhou HP, Ma CL, Wang PP. Regulation of plant responses to salt stress. Int J Mol Sci. 2021;22(9):4609.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Ma LY, Zhang H, Sun LR, Jiao YH, Zhang GZ, Miao C, et al. NADPH oxidase AtrbohD and AtrbohF function in ROS-dependent regulation of Na+/K+ homeostasis in Arabidopsis under salt stress. J Exp Bot. 2012;63(1):305–17.

    Article  CAS  PubMed  Google Scholar 

  36. Yang YL, Xu SJ, An LZ, Chen NL. NADPH oxidase-dependent hydrogen peroxide production, induced by salinity stress, may be involved in the regulation of total calcium in roots of wheat. J Plant Physiol. 2007;164(11):1429–35.

    Article  CAS  PubMed  Google Scholar 

  37. Cao BL, Li N, Xu K. Crosstalk of phenylpropanoid biosynthesis with hormone signaling in Chinese cabbage is key to counteracting salt stress. Environ Exp Bot. 2020;179:104209.

    Article  CAS  Google Scholar 

  38. Goujon T, Sibout R, Pollet B, Maba B, Nussaume L, Bechtold N, et al. A new Arabidopsis thaliana mutant deficient in the expression of O-methyltransferase impacts lignins and sinapoyl esters. Plant Mol Biol. 2003;51(6):973–89.

    Article  CAS  PubMed  Google Scholar 

  39. Sewalt VJH, Ni WT, Jung HG, Dixon RA. Lignin impact on fiber degradation: increased enzymatic digestibility of genetically engineered tobacco (Nicotiana tabacum) stems reduced in lignin content. J Agric Food Chem. 1997;45(5):1977–83.

    Article  CAS  Google Scholar 

  40. Zhang K, Cui HT, Cao SH, Yan L, Li MN, Sun Y. Overexpression of CrCOMT from Carex rigescens increases salt stress and modulates melatonin synthesis in Arabidopsis thaliana. Plant Cell Rep. 2019;38(12):1501–14.

    Article  CAS  PubMed  Google Scholar 

  41. Gurmani AR, Bano A, Khan SU, Din J, Zhang JL. Alleviation of salt stress by seed treatment with abscisic acid (ABA), 6-benzylaminopurine (BA) and chlormequat chloride (CCC) optimizes ion and organic matter accumulation and increases yield of rice (Oryza sativa L.). AJCS. 2011;5(10):1278–85.

    CAS  Google Scholar 

  42. Afzal I, Basra SMA, Farooq M, Nawaz A. Alleviation of salinity stress in spring wheat by hormonal priming with ABA, salicylic acid and ascorbic acid. Int J Agri Biol. 2006;8(1):23–8.

    CAS  Google Scholar 

  43. Gurmani AR, Bano A, Salim M. Effect of abscisic acid and benzyladenine on growth and ion accumulation of wheat under salinity stress. Pak J Bot. 2007;39(1):141–9.

    Google Scholar 

  44. Saeedipour S. Salinity tolerance of rice lines related to endogenous abscisic acid (ABA) level synthesis under stress. African J Plant Sci. 2011;5(11):628–33.

    CAS  Google Scholar 

  45. Moons A, Keyser AD, Montagu MV. A group 3 LEA cDNA of rice, responsive to abscisic acid, but not to jasmonic acid, shows variety-specific differences in salt stress response. Gene. 1997;191(2):197–204.

    Article  CAS  PubMed  Google Scholar 

  46. Lee SK, Sohn EY, Hamayun M, Yoon JY, Lee IJ. Effect of silicon on growth and salinity stress of soybean plant grown under hydroponic system. Agrofor Syst. 2010;80(3):333–40.

    Article  Google Scholar 

  47. Mambro RD, Ruvo MD, Pacifici E, Salvi E, Sozzani R, Benfey PN, et al. Auxin minimum triggers the developmental switch from cell division to cell differentiation in the Arabidopsis root. PNAS. 2017;114(36):7641–9.

    Article  Google Scholar 

  48. Garrido-Vargas F, Godoy T, Tejos R, O'Brien JA. Overexpression of the auxin receptor AFB3 in Arabidopsis results in salt stress resistance and the modulation of NAC4 and SZF1. Int J Mol Sci. 2020;21(24):9528.

    Article  CAS  PubMed Central  Google Scholar 

  49. Albert B, Kos-Braun IC, Henras AK, Dez C, Rueda MP, Zhang X, et al. A ribosome assembly stress response regulates transcription to maintain proteome homeostasis. eLife. 2019;8:e45002.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Rippmann JF, Michalowski CB, Nelson DE, Bohnert HJ. Induction of a ribosome-inactivating protein upon environmental stress. Plant Mol Biol. 1997;35:701–9.

    Article  CAS  PubMed  Google Scholar 

  51. Wu CCC, Peterson A, Zinshteyn B, Regot S, Green R. Ribosome collisions trigger general stress responses to regulate cell fate. Cell. 2020;182(2):404–16.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Iordanov MS, Pribnow D, Magun JL, Dinh TH, Pearson JA, Chen SL, et al. Ribotoxic stress response: activation of the stress-activated protein kinase JNK1 by inhibitors of the peptidyl transferase reaction and by sequence-specific RNA damage to the alpha-sarcin/ricin loop in the 28S rRNA. Mol Cell Biol. 1997;17(6):3373–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30(15):2114–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Langmead B, Wilks C, Antonescu V, Charles R. Scaling read aligners to hundreds of threads on general-purpose processors. Bioinformatics. 2019;35(3):421–32.

    Article  CAS  PubMed  Google Scholar 

  55. Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14(4):R36.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550.

    Article  PubMed  PubMed Central  Google Scholar 

  57. Leng N, Dawson JA, Thomson JA, Ruotti V, Rissman AI, Smits BMG, et al. EBSeq: an empirical Bayes hierarchical model for inference in RNA-seq experiments. Bioinformatics. 2013;29(8):1035–43.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Ju WY, Li JX, Yu WR, Zhang RC. iGraph: an incremental data processing system for dynamic graph. Front Comput Sci. 2016;10:462–76.

    Article  Google Scholar 

  59. Morris JH, Apeltsin L, Newman AM, Baumbach J, Wittkop T, Su G, et al. clusterMaker: a multi-algorithm clustering plugin for Cytoscape. BMC Bioinformatics. 2011;12:436.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Kanehisa M, Goto S, Kawashima S, Nakaya A. The KEGG databases at GenomeNet. Nucleic Acids Res. 2002;30(1):42–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Wei W, Tao JJ, Chen HW, Li QT, Zhang WK, Ma B. A histone code reader and a transcriptional activator interact to regulate genes for salt tolerance. Plant Physiol. 2017;175(3):1304–20.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Lisso J, Schröder F, Fisahn J, Müssig C. NFX1-LIKE2 (NFXL2) suppresses abscisic acid accumulation and stomatal closure in Arabidopsis thaliana. Plos One. 2011;6(11):e26982.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. Zhang XX, Liu SK, Takano T. Two cysteine proteinase inhibitors from Arabidopsis thaliana, AtCYSa and AtCYSb, increasing the salt, drought, oxidation and cold tolerance. Plant Mol Biol. 2008;68(1–2):131–43.

    Article  CAS  PubMed  Google Scholar 

  64. Liu J, Pang X, Cheng Y, Yin YH, Zhang Q, Su WB, et al. The Hsp70 gene family in Solanum tuberosum: genome-wide identification, phylogeny, and expression patterns. Sci Rep. 2018;8(1):16628.

    Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank Ran Xu (Shandong Academy of Agricultural Sciences) and Ying Wang (Jilin University) for the materials of QH34 and DN50, respectively.

Statement

All methods were carried out in accordance with relevant guidelines and regulations.

Funding

This work was supported by the National Key Research and Development Program (grant no. 2021YFF1001100), the Taishan Scholars Program of Shandong Province (tsqn201812036), the Agricultural Variety Improvement Project of Shandong Province (2019LZGC004), the Program for Scientific Research Innovation Team of Young Scholar in Colleges and Universities of Shandong Province, China (2020KJF008) and Natural Science Foundation of Shandong Province (ZR2020QC219).

Author information

Author notes

  1. Junmei Hu and Yongbin Zhuang contributed equally to this work.

Authors and Affiliations

  1. College of Agronomy, State Key Laboratory of Crop Biology, Shandong Agricultural University, Tai’an, 271018, Shandong, China

    Junmei Hu, Yongbin Zhuang, Xianchong Li, Xiaoming Li & Dajian Zhang

  2. Crop Research Institute, Shandong Academy of Agricultural Sciences, Ji’nan, 250131, Shandong, China

    Ran Xu

  3. College of Life Sciences, Yantai University, Yan’tai, 264005, Shandong, China

    Chanchan Sun

  4. The Key Laboratory of Plant Development and Environmental Adaptation Biology, Ministry of Education, School of Life Sciences, Shandong University, Qing’dao, 266237, Shandong, China

    Zhaojun Ding

Authors
  1. Junmei Hu
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Yongbin Zhuang
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Xianchong Li
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Xiaoming Li
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Chanchan Sun
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Zhaojun Ding
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Ran Xu
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Dajian Zhang
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

RX and DJZ conceived and designed the research. JMH, XCL, XML, CCS and ZJD performed the experiment. YBZ analyzed the transcriptome data. YBZ and JMH drafted the manuscript. All authors read and approved the final version of the manuscript.

Corresponding authors

Correspondence to Ran Xu or Dajian Zhang.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Fig. S1.

The measurement of salt tolerance in 22 soybean accessions.

Additional file 2: Fig. S2.

Evaluation of RNA-seq data from the two soybean cultivars DN50 and QH34.

Additional file 3: Fig. S3.

Heatmap of negative regulation abscisic acid activity signaling pathway.

Additional file 4: Fig. S4.

Heatmap of oxidoreductase activity signaling pathway.

Additional file 5: Fig. S5.

Heatmap of negative regulation of proteolysis signaling pathway.

Additional file 6: Fig. S6.

GO enrichment analysis for uDEGs.

Additional file 7: Fig. S7.

KEGG pathway analysis of uDEGs.

Additional file 8: Fig. S8.

Pathway analysis of gene network clusters for sDEGs.

Additional file 9: Fig. S9.

Heatmap of hormone signaling pathways.

Additional file 10: Fig. S10.

Heatmap of MAPK signaling pathway.

Additional file 11: Fig. S11.

Heatmap of phenylpropanoid biosynthesis signaling pathway.

Additional file 12: Fig. S12.

Heatmap of ribosome metabolism biosynthesis signaling pathway.

Additional file 13: Table S1.

The source information of 22 soybean germplasms.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Hu, J., Zhuang, Y., Li, X. et al. Time-series transcriptome comparison reveals the gene regulation network under salt stress in soybean (Glycine max) roots. BMC Plant Biol 22, 157 (2022). https://doi.org/10.1186/s12870-022-03541-9

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12870-022-03541-9

Keywords

BMC Plant Biology

ISSN: 1471-2229

Contact us