Integrated transcriptome and proteome analysis reveals complex regulatory mechanism of cotton in response to salt stress

Abstract

Background

Soil salt stress seriously restricts the yield and quality of cotton worldwide. To investigate the molecular mechanism of cotton response to salt stress, a main cultivated variety Gossypium hirsutum L. acc. Xinluzhong 54 was used to perform transcriptome and proteome integrated analysis.

Results

Through transcriptome analysis in cotton leaves under salt stress for 0 h (T0), 3 h (T3) and 12 h (T12), we identified 8 436, 11 628 and 6 311 differentially expressed genes (DEGs) in T3 vs. T0, T12 vs. T0 and T12 vs. T3, respectively. A total of 459 differentially expressed proteins (DEPs) were identified by proteomic analysis, of which 273, 99 and 260 DEPs were identified in T3 vs. T0, T12 vs. T0 and T12 vs. T3, respectively. Metabolic pathways, biosynthesis of secondary metabolites, photosynthesis and plant hormone signal transduction were enriched among the identified DEGs or DEPs. Detail analysis of the DEGs or DEPs revealed that complex signaling pathways, such as abscisic acid (ABA) and jasmonic acid (JA) signaling, calcium signaling, mitogen-activated protein kinase (MAPK) signaling cascade, transcription factors, activation of antioxidant and ion transporters, were participated in regulating salt response in cotton.

Conclusions

Our research not only contributed to understand the mechanism of cotton response to salt stress, but also identified nine candidate genes, which might be useful for molecular breeding to improve salt-tolerance in cotton.

Background

Soil salinity caused by high concentration of salt ions in soil, which seriously restricts the agricultural production. More than 800 million hectares of land and 20% of the arable land throughout the world are affected by salt stress, which will continue to deteriorate (Al Murad et al. 2020; Munns and Tester 2008). It was showed that the problem of soil salinity become more and more serious in China. The area of saline-alkali land in China is about 100 million hectares. Moreover, the area of secondary salinization in the western region and the Yellow River Delta are increased 150 000∼200 000 and 100 000 hm² every year, respectively. The high concentration of salt ions can change soil structure, reduce soil permeability and water conductivity, trigger different levels of stress response on plants, such as osmotic stress, ion cytotoxicity, nutritional deficiency, and oxidative stress (Mahajan and Tuteja 2005; Muchate et al. 2016; Zelm et al. 2020). Osmotic stress reduces water uptake and leads to water deficit in plants, which occurs in the early stage of salt stress. Ion cytotoxicity is caused by the high concentration of salt ions (Na⁺) in the cytoplasm, and it disrupts the uptake of other ions and obstructs some metabolic pathways, such as photosynthesis (Wang et al. 2020; Zhang et al. 2019). Osmotic stress and ionic stress can cause secondary stress in plants, for example, the high concentration of reactive oxygen species will cause damage to cell structures and biological macromolecules (Yang and Guo 2018a).

The salt tolerance of plants depends on their ability to extract water and nutrients from saline soils and to avoid excessive accumulation of salt ions in plant tissues (Zorb et al. 2019). Most of crop plants are salt-sensitive, and the yield of crops will be seriously reduced when they are subjected to moderate salt stress (Park et al. 2016; Zorb et al. 2019). Therefore, analyzing the mechanism of plant response to salt stress and cultivating new salt-tolerant crop varieties are very important to ensure the sustainability of agricultural production and food security.

Cotton (Gossypium) is the most important fiber crop, and China is a larger producer and the largest consumer of cotton in the world. Xinjiang is the main cotton-producing area in China, and more than 32.6% of cultivated land in Xinjiang is affected by different degrees of salinization stress, which seriously threatens the safety of cotton production. It is very important to improve salt tolerance based on ensuring high yield and quality in cotton. Transcriptome and proteome have made progress in revealing the mechanism of salt tolerance and identifying candidate genes in cotton (Gong et al. 2017; Guo et al. 2015; Li et al. 2015; Peng et al. 2014; Shan et al. 2019; Sikder et al. 2020; Yang et al. 2019). As there are multi-level regulatory machineries exist in salt-stress response, including transcription and translation regulations, it is important to monitor the gene expression level of RNA and protein simultaneously. Fortunately, the development of integrated transcriptome and proteome makes this research strategy possible (Chen et al. 2015; Peng et al. 2018; Trevisan et al. 2015; Wang et al. 2014; Yang et al. 2019).

In the present study, an integrated transcriptome and proteome analysis were performed to deeply investigate the molecular mechanism of cotton in response to salt stress. There were 15 822 DEGs and 459 DEPs identified in mRNA transcription and protein expression level, respectively, and a total of 164 stress-response associated DEGs and DEPs were identified. And important biological pathways related to salt-stress response were elucidated. Furthermore, 9 candidate genes were identified by integrating proteomic and transcriptomic profiles. Our results will further enrich the understanding of molecular mechanism in cotton in response to salt stress.

Materials and methods

Plant materials and treatments

The widely planted upland cotton variety Xinluzhong 54 (Gossypium hirsutum L. acc.) was used in the research, and it was provided by Xinjiang Academy of Agricultural Sciences (Additional file 1: Figure S1). Cotton seedlings were cultured in Hoagland solution under 16 h light/8 h dark conditions. For salt stress treatment, 250 mmol·L^− 1 NaCl was added to the solution at two-leaf stage seedlings to simulate salt stress. The second newly developed leaves were collected at 0 h (T0), 3 h (T3) and 12 h (T12) for subsequent RNA and protein isolations.

RNA sequencing (RNA-seq) and isobaric tags for relative and absolute quantification (iTRAQ) analysis

Total RNAs were extracted using the RNAprep Pure Plant Kit (Cat. #DP441, TIANGEN). RNA-Seq was conducted using Illumina Hiseq platform (BGI Tech, Shenzhen, China). After removing the low-quality, adaptor-polluted and high content of unknown base (N) reads, the filtered reads were mapped against Gossypium hirsutum (acc. TM-1) genome and genes were annotated accordingly (Zhang et al. 2015). The gene expression levels for each sample were calculated with RSEM software, the DEGs were detected with DEseq2 (Fold change ≥2, adjusted P-value ≤0.05).

For iTRAQ analysis, the protein quantitative analysis is carried out by BGI company (BGI Tech, Shenzhen, China), and the main processes includes: protein extraction and digestion, peptide labeling and fractionation, protein quantification by LC-MS / MS. The quantitative analysis was performed by IQuant software (Wen et al. 2014). T3 vs. T0 (T3/T0), T12 vs. T0 (T12/T0) and T12 vs. T3 (T12/T3) were set as comparison groups, the DEPs were identified using the threthold of fold change > 1.2 (or < 0.83) and Q-value < 0.05. Both RNA-Seq and iTRAQ experiments were performed with three biological replicates.

Functional enrichment analysis

KEGG (Kyoto encyclopedia of genes and genomes) and GO (Gene ontology) functional enrichment analysis of DEGs and DEPs were performed using the “Gene-list enrichment” program in KOBAS3.0 database (http://kobas.cbi.pku.edu.cn/anno_iden.php) (Xie et al. 2011). Transcription factors were identified by PlantTFDB 4.0 (http://planttfdb.cbi.pku.edu.cn/). The protein interaction networks were predicted by STRING software (ver. 11.0), and were visualized by Cytoscape software (ver. 3.4.0). The gene/protein expression patterns were showed by heatmap with the expression values normalized by Genesis software.

Quantitative real-time PCR (qRT-PCR) analysis

High-quality cDNA was obtained by reverse transcription using SuperScript III Reverse Transcriptase following the manufacturer’s instructions (Cat. No.18080–093, Invitrogen). Gene expression levels were detected using an ABI Prism 7500 system (Applied Biosystems). GhUBQ7 (GenBank accession No. DQ116441) was used as the internal control. Gene-specific primers for qRT-PCR were designed according to the cDNA sequences using Primer Premier 5.0 software and synthesized commercially (Genscript Bioscience) (Additional file 1: Table S1).

Results

Identification of DEGs and DEPs during cotton response to salt stress

The seedlings of Xinluzhong 54 (Gossypium hirsutum L. acc.) were cultured in Hoagland solution and 250 mmol·L^− 1 NaCl was added to simulate salt stress. At the two-leaf stage, the second newly-developed leaves were used for transcriptome and proteome analysis, the workflow of experimental design and analysis was shown in Fig. 1a. Three samples with three biological replicates were analyzed using Illumina Hiseq platform and iTRAQ system.

Fig. 1

Identification of differentially expressed genes (DEGs) and proteins (DEPs) in response to salt stress. a Scheme of experimental design in integrated transcriptome and proteome analysis. b Number of DEGs in different comparison groups. c Venn diagram of DEGs. d Number of DEPs in different comparison groups. e Venn diagram of DEPs

The average output data of each sample from RNA-Seq was 6.64 Gb. A Phred quality score (Q score) was used to evaluate the base call accuracy of the read data. Q20 represents an incorrect base call of 1 in 100 (99%), Q30 represents an incorrect base call of 1 in 1000 (99.9%) and Q30 is usually used as an evaluation standard in the quality control of sequencing data. Here, the average Q20 and Q30 of each sample were 97.77 and 93.65%, respectively. The clean reads were mapped to the reference genome of Gossypium hirsutum (TM-1), and the average alignment rate of each sample is 82.83% (Additional file 1: Table S2). DEseq2 algorithm was used to detect the DEGs with the standard of fold change ≥2 and P-value ≤0.05. There were 8 436 (4 123 up-regulated, 4 313 down-regulated), 11 628 (4 649 up-regulated, 6 979 down-regulated) and 6 311 (2 169 up-regulated, 4 142 down-regulated) DEGs were identified in T3/T0, T12/T0 and T12/T3, respectively (Fig. 1b). DEGs were distributed more in D sub-genome (53%) than A sub-genome. And more DEGs were distributed on chromosomes A05 and D05, while less were distributed on A04 and D03 (Additional file 1: Figure S2a). Twenty-five genes were randomly selected to evaluate the RNA-Seq results by qRT-PCR, there was a significant correlation between these two groups of data (R² = 0.874) (Additional file 1: Figure S2b). It indicated that the sequencing data of the RNA-Seq had high reliability, and could be used for subsequent analysis. Venn diagram showed that 5 402, 3 945 and 2 342 common DEGs were detected in T3 / T0 and T12 / T0, T12 / T0 and T12 / T3, T3 / T0 and T12 / T3, respectively, with 1 136 common DEGs shared by the three comparisons (Fig. 1c). The 1 136 common DEGs showed up-regulated and down-regulated expression patterns by calculating hierarchical clustering analysis, with up-regulated genes enriching in galactose metabolism, biosynthesis of secondary metabolites, metabolic pathways, valine, leucine and isoleucine degradation, glucosinolate biosynthesis, and down-regulated genes in metabolic pathways, carbon fixation in photosynthetic organisms, biosynthesis of secondary metabolites, carbon metabolism and flavonoid biosynthesis (Additional file 1: Figure S3a).

A total of 459 DEPs were identified by iTRAQ, and 273 (185 up-regulated, 88 down-regulated), 99 (38 up-regulated, 61 down-regulated) and 260 (106 up-regulated, 154 down-regulated) DEPs were identified in T3/T0, T12/T0 and T12/T3, respectively (Fig. 1d). All DEPs encoding genes were equally distributed on A sub-genome and D sub-genome, and more DEPs were distributed on chromosomes A05 and D05. Venn diagram showed that 18, 127 and 35 common DEPs in T3/T0 and T12/T0, T3/T0 and T12/T3, T12/T0 and T12/T3, respectively (Fig. 1e). There were seven common DEPs in three groups, three proteins (Gh_A02G0551, Gh_A07G0867, Gh_D08G1902) were up-regulated at 3 h, two proteins (Gh_D02G2126, Gh_A03G1705) were down-regulated at 3 h and up-regulated at 12 h, and Gh_Sca010764G01 was up-regulated at 3 h and 12 h, Gh_D11G1672 was down-regulated at 3 h and 12 h.

Integration of proteomic and transcriptomic profiles

To investigate the multi-level regulation of gene expression under salt stress, an integrative analysis of the transcriptome and proteome was performed in this study. Results showed that 77, 70 and 47 associated DEGs and DEPs were identified in T3/T0, T12/T0 and T12/T3, respectively (Fig. 2a, Additional file 1: Table S3). In T3/T0, 31 and 46 genes showed the same and opposite expression patterns at the transcriptional level and translational level, respectively. In T12/T0, 62 and 8 genes showed the same and opposite expression patterns at these two levels, respectively. In T12/T3, 31 and 16 genes showed the same and opposite expression patterns at these two levels, respectively (Fig. 2a). These results suggested that genes in response to salt stress may undergo complex regulation at the transcriptional and translational levels.

Fig. 2

Statistics and functional analysis on the correlation between DEGs and DEPs. a Statistics on the number of correlated DEGs and DEPs. b Venn diagram of correlated DEGs/DEPs. c Concordance between the gene changes at transcriptome and proteome levels. d Functional enrichment analysis of the correlated DEGs and DEPs

Venn diagram showed that nine common, associated genes (Gh_D08G1902, Gh_A07G0867, Gh_Sca010764G01, Gh_D13G1729, Gh_D09G1072, Gh_D06G2351, Gh_A11G0835, Gh_D02G2126, Gh_A03G1705) were found in the comparion group of T3/T0 and T12/T0. Eleven associated genes (Gh_D11G1672, Gh_D10G2299, Gh_D06G1799, Gh_D05G2157, Gh_A11G1830, Gh_A10G1991, Gh_A06G1466, Gh_A05G1452, Gh_A01G1839, Gh_D02G2126, Gh_A03G1705) were identified in T12/T0 and T12/T3. And twelve associated genes (Gh_D09G1409, Gh_D06G0164, Gh_D04G0628, Gh_D02G0438, Gh_A04G0778, Gh_A03G1240, Gh_A01G1605, Gh_D05G0122, Gh_A05G1494, Gh_A03G0173, Gh_D02G2126, Gh_A03G1705) were identified in T3/T0 and T12/T3 (Fig. 2b). All associated genes were listed in Additional file 1: Table S4. Spearman correlation analysis showed that a poor correlation between the expression changes of these associated DEGs and DEPs was observed (R: − 0.37∼ – 0.71) (Fig. 2c).

Functional enrichment analysis of DEGs and DEPs

KEGG analysis was performed to study the functional enrichment of up-regulated and down-regulated DEGs in each comparison group. For up-regulated DEGs, 21, 30 and 4 pathways were enriched in T3/T0, T12/T0 and T12/T3, respectively (Additional file 1: Figure S3b). Plant hormone signal transduction, biosynthesis of secondary metabolites and metabolic pathways occurred simultaneously in these three groups. For down-regulated DEGs, 17, 22 and 25 pathways were enriched in T3/T0, T12/T0 and T12/T3, respectively (Additional file 1: Figure S3c). Ribosome was the most significantly enriched pathway in T3/T0, and metabolic pathways was significantly enriched in T12/T0 and T12/T3.

For up-regulated DEPs, 8, 2 and 14 pathways were enriched in T3/T0, T12/T0 and T12/T3, respectively (Corrected P ≤ 0.01) (Additional file 1: Figure S4a). Photosynthesis and pyruvate metabolism were the common pathways in these three groups. Glyoxylate and dicarboxylate metabolism, ribosome, carbon fixation in photosynthetic organisms and citrate cycle (TCA cycle) were significantly enriched in T3/T0. In addition, ten pathways were enriched in T12/T3, such as endocytosis, alpha-linolenic acid metabolism and glycerophospholipid metabolism. For the down-regulated DEPs, 11, 13 and 9 pathways were enriched in T3/T0, T12/T0 and T12/T3, respectively (Additional file 1: Figure S4b). Metabolic pathways and biosynthesis of secondary metabolites were the common pathways in these three groups. In addition, selenocompound metabolism, biosynthesis of amino acids and porphyrin and chlorophyll metabolism were commonly enriched in T3/T0 and T12/T0, photosynthesis and ribosome were commonly enriched in T3/T0 and T12/T3, flavonoid biosynthesis and carbon fixation in photosynthetic organisms were commonly enriched in T12/T0 and T12/T3.

KEGG analysis showed that these associated genes/proteins were enriched in multiple pathways (Corrected P-value ≤0.05) (Fig. 2d). We found some DEGs showed the same changing patterns at transcriptional and protein levels and were enriched in some pathways, such as flavonoid biosynthesis (Gh_D08G1902, Gh_A10G1079, Gh_A05G1424, Gh_A12G0549), porphyrin and chlorophyll metabolism (Gh_A09G1143, Gh_A05G3766, Gh_D01G1828, Gh_A10G0282, Gh_A10G2274), thiamine metabolism (Gh_A06G1269, Gh_A05G1432, Gh_D05G1605), biosynthesis of amino acid (Gh_D08G2317, Gh_A13G1464, Gh_D13G1473, Gh_D06G1578, Gh_D01G1099, Gh_D05G2157, Gh_D09G0576).

Disturbance of energy metabolism caused by salt stress

The carbohydrate metabolism pathways and photosynthesis were very important for the storage and utilization of energy in plants. Here, the changes of DEGs and DEPs which were involved in these two biological processes were analyzed. Treated with salt stress, cotton genes involved in the glycolysis/ gluconeogenesis pathways exhibited strong expression levels. The expression levels of GhPFK3 (Gh_A05G0198), GhPKP2 (Gh_A08G0871) and GhPKP1 (Gh_A10G1036) were up-regulated under salt stress. GhHXK3 (Gh_D06G0003), a gene that encodes hexokinase (HK), was down-regulated under salt stress, with the transcription levels decreased to 22 and 9.6% of control at 3 h and 12 h, respectively (Fig. 3a). The protein abundance of enolase PGH1 (Gh_Sca286293G01), fructose-bisphosphate aldolase FBA5 (Gh_D13G1361), pyruvate dehydrogenase subunit (Gh_A12G1299) and triosephosphate isomerase (Gh_D01G1600) increased significantly at 3 h after salt stress (Additional file 1: Figure S5a). The acetaldehyde dehydrogenase protein (Gh_D06G1578) increased significantly at 12 h, hexokinase HK (Gh_A13G1742) and ethanol dehydrogenase (Gh_A01G1605) protein decreased at 3 h.

Fig. 3

The expression patterns of functional genes which were involved in energy metabolism. a The expression patterns of genes involved in the glycolysis/gluconeogenesis pathways. b The expression patterns of genes involved in tricarboxylic acid cycle. c The expression patterns of genes involved in photosynthesis

The down-regulation was observed in DEGs encoding isocitrate dehydrogenase (IDH) (GhCICDH:, Gh_A11G1562), pyruvate dehydrogenase E1 subunit (GhMAB1:, Gh_A12G1299) and malate dehydrogenase (GhMMDH2: Gh_A01G04040, GhMMDH1: Gh_A04G0320) which were involved in the tricarboxylic acid (TCA) cycle (Fig. 3b). And the abundance of malate dehydrogenase (Gh_D02G0438, Gh_A02G0386, Gh_D05G3328) and pyruvate dehydrogenase E1 subunit (Gh_A12G1299) increased at 3 h after salt stress (Additional file 1: Figure S5b). The protein abundance of succinate dehydrogenase subunit (Gh_A10G1110) and 2-oxoglutarate dehydrogenase (Gh_A05G3057) decreased at 3 h (Additional file 1: Figure S5b).

Similarly, the expression levels of genes modulating photosynthesis was significantly downregulated in cotton. As photosystem I subunit PSAO related gene, GhPSAO (Gh_D07G1090) and light-harvesting chlorophyll protein complex LHC-related genes, GhLHCB1 (Gh_D06G2351), GhLHCB2 (Gh_A07G1725) and GhLHCA4 (Gh_A12G1617) (Fig. 3c). The protein abundances of most DEPs involved in photosynthesis were increased under salt stress (Additional file 1: Figure S5c). For example, photosystem II subunit PSBO (Gh_D11G1897), PSBP (Gh_A05G3293, Gh_D04G0311) and PSBQ2 (Gh_A10G2195) proteins increased at 3 h under salt stress. The abundance of photosystem I related proteins psaA (Gh_A01G1364, Gh_A09G0964) and psaB (Gh_A12G1136, Gh_Sca005646G01) decreased after 3 h of salt stress, but increased at 12 h. Taken together, these results suggested that salt stress can cause the changes in expression level of cotton genes involved in carbohydrate metabolism pathways and photosynthesis pathways, to enhance the tolerance and survival of cotton in salt stress.

Signal transduction in response to salt stress in cotton

Plant hormone signal transduction

Plant hormones play an important role in the regulation of plant salt stress resistance. Here, the expression changes of key genes in different hormone signaling pathways were analyzed (Figure S6). ABA and JA hormone signals are important for plants to resist salt stress. Most of the genes involved in ABA signaling pathway, such as PP2C, SnRK2 and ABF, were up-regulated by salt stress (Fig. 4a). We speculated that the continuous activation of ABA signal may play an important role in cotton response to salt stress. ABF transcription factors play an important role in regulating ABA signal transduction. Fourteen genes encoding ABF transcription factors were identified, all of them were up-regulated by salt stress (Fig. 4a). Nine ABF genes (Gh_D05G2495, Gh_A05G2234, Gh_D12G0214, Gh_A05G1751, Gh_D05G1946, Gh_D02G0909, Gh_A03G1664, Gh_A03G2095, Gh_D02G2079) were continuously induced by salt stress during 0–12 h, and the maximum expression levels of four genes (Gh_A12G0212, Gh_D13G0400, Gh_A13G0355, Gh_A13G2006) were detected at 3 h (Fig. 4a).

Fig. 4

Plant signal transduction in response to salt stress in cotton. a-b The expression patterns of genes which were involved in ABA and JA signal pathways. c The expression patterns of CBL and CIPK genes. d The expression patterns of MAPK signaling cascade genes

The JAR1 gene in JA pathway encodes a jasmonic acid-amino synthase, which could catalyze the formation of bioactive jasmonic acid-isoleucine (JA-Ile) in Arabidopsis. The JAR1 gene in cotton, Gh_A08G1120 and Gh_D08G1403, were significantly up-regulated in T3/T0, so we speculated that salt stress might induce JA-Ile synthesis and activate JA signal in the early stage (Fig. 4b). MYC2 is an important transcription factor in JA signaling pathway in Arabidopsis. We identified two genes encoding MYC2 transcription factors that were differentially expressed, Gh_D08G1707 was significantly induced by salt stress and showed continuous expression patterns, while Gh_D09G1895 was down-regulated by salt stress. Furthermore, 21 differentially expressed JAZ genes were identified, among them, 18 genes were up-regulated (Fig. 4b).

Transcription factor

Transcription factor plays an important role in plant response to abiotic stress. There were 696 (up-regulated: 429, down-regulated: 267), 926 (up-regulated: 492, down-regulated: 434) and 532 (up-regulated: 232, down-regulated: 300) transcription factors identified in T3/T0, T12/T0 and T12/T3, respectively (Figure S7a and b). Among these transcription factors, NAC, ERF, MYB, bZIP and WRKY family genes accounted for a large proportion. Among the transcription factors with the spatio-temporal expression patterns, about 32 and 49% were only differentially expressed at 3 h and 12 h after salt stress, respectively (Additional file 1: Figure S7a and b). Some of the transcription factors that were revealed to be responsive to salt stress previous studies were also detected in this study (Gao et al. 2020; Krishnamurthy et al. 2020; Xie et al. 2010; Zhao et al. 2019) (Additional file 1: Figure S7c).

Calcineurin B-like (CBL)-interacting protein kinase (CIPK) signaling cascade

The CBL-CIPK signaling system can participate in the process of plant stress response by sensing calcium ion to regulate intracellular signal transduction. Two genes encoded CBL, GhCBL2 (Gh_A03G0846) and GhCBL3 (Gh_A13G1099), were identified, both were significantly up-regulated by salt stress. There were 43 CIPK genes identified, and most of these genes were up-regulated by salt stress (Fig. 4c). We predicted the potential protein interaction models between these CBL and CIPK proteins through the STRING database. Results showed that these two CBL proteins could interact with multiple CIPK proteins, respectively (Additional file 1: Figure S8a). Furthermore, some CBL and CIPK genes were found with similar expression patterns under salt stress, such as CBL2 and CIPK8 (Gh_A05G3756, Gh_A06G1872, Gh_D05G2084), both types were continuously induced by salt stress.

Mitogen-activated protein kinase signaling cascade

Mitogen-activated protein kinase (MAPK) cascade signaling pathway is composed of three kinds of serine/threonine protein kinases, MAP3K, MKK and MPK. There were 18, 4 and 17 genes encoded MAP3K, MKK, MPK identified, respectively (Fig. 4d). Most of the MAPK pathway genes were up-regulated by salt stress at 3 h. The potential protein interaction models were also predicted through the STRING database. Two up-regulated genes, Gh_A07G0124 and Gh_D07G2384, which encoding GhMKK2, were predicted wide interaction models with upstream MAP3K proteins and downstream MPK proteins. We speculated that GhMKK2 might play an important role in regulating cotton response to salt stress through MAPK cascade pathway (Additional file 1: Figure S8b).

Antioxidant and ion transport systems

Salt stress causes the imbalance of reactive oxygen species (ROS) homeostasis and activate the ROS scavenging mechanism. NADPH oxidase is involved in the production of ROS in plants under stress. Four genes (GhRBOHA; Gh_A02G1791, Gh_A12G2669, Gh_D03G0688 and Gh_D12G2750) encoding NADPH oxidase respiratory burst oxidase homolog A (RBOHA) were up-regulated by salt stress (Fig. 5a). Consistently, dramatic up-regulation was detected in genes encoding catalase (GhCAT2), superoxide dismutase (GhFSD2) and glutathione S-transferase (GhGSTT1, GhGSTZ1) (Fig. 5b). Our results indicated that pathways involved of CAT, SOD and GST genes might play an important role in protecting cotton from oxidative damage caused by salt stress.

Fig. 5

The expression patterns of genes involved in antioxidant and ion transport systems. a-b Expression patterns of genes related to reactive oxygen species production and scavenging. c-d Expression patterns of ion transport related genes

The regulation of intracellular ion homeostasis is very important for plants response to salt stress. Ion channel protein (Na⁺/H⁺/K⁺ transporter) are tightly linked to plant response to salt stress. Thus, we analyzed the expression of genes that encoding ion transporters. GhNHX2 (Gh_A11G2132) and GhKUP3 (Gh_D04G0700), GhKUP11 (Gh_A12G0438) encoding Na⁺/H⁺ transporters and K⁺ transporters, were significantly up-regulated under salt stress (Fig. 5c). GhHKT1 (Gh_A08G2473) encoding Na⁺ transporters were down-regulated under salt stress at 12 h. Potassium channels are involved in the regulation of ion absorption and transport. And significant up-regulation at 12 h under salt stress were observed in genes that encoding inward rectifying potassium channel (GhAKT1, Gh_D13G2111) and outward rectifying K⁺ channel (GhTPK1, Gh_D01G0726; GhKCO5, Gh_A11G0947; GhSKOR, Gh_D05G1041) (Fig. 5d).

The identification of candidate genes by integrating proteomic and transcriptomic profiles

Based on the gene expression and functional annotation, we identified 9 candidate genes that might play an important role in cotton response to salt stress, i.e., GhLEA14 (Gh_A11G0835), GhKCS6 (Gh_A03G1286), GhPRXR1 (Gh_A05G1452), GhPOD (Gh_A08G0714), GhFIB (Gh_A05G1494), GhALDH7B4 (Gh_D06G1578), GhUSP (Gh_A03G0386) and GhAVP1 (Gh_D05G0122), GhADH1 (Gh_A01G1605) (Table 1). GhLEA14 encodes late embryogenesis abundant protein and its expression level was up-regulated at both transcriptional and protein levels under salt stress. This gene might be involved in the regulation of osmotic stress and the prevention of oxidant stress in cotton undersalt stress. The peroxidase superfamily protein GhPRXR1 and GhPOD were up-regulated at transcription and protein levels by salt stress, two correspondent genes might be involved in the hydrogen peroxide catabolic process. GhADH1 and GhALDH7B4 encoding ethanol dehydrogenase and acetaldehyde dehydrogenase, respectively. These two genes were up-regulated at the transcriptional and protein levels under salt stress. We speculated that these two genes may be involved in the cellular detoxification of cotton under salt stress.

Table 1 Candidate genes involved in salt stress response in cotton

Discussion

Complex regulatory mechanisms took place when plants encountered salt stress. Although cotton is an inherently salt-tolerant crop, with the aggravation of soil salinization in its planting area, it is important to analyze the molecular mechanism of its response to salt stress and to cultivate novel cotton varieties with higher salt tolerance. Due to the high yield, perfect fiber quality and disease resistance, Xinluzhong 54 is one of the main cotton variety grown in Xinjiang nowadays (Additional file 1: Figure S1a). However, soil salinity is one of the major factors restricting the further expansion of Xinluzhong 54 (Additional file 1: Figure S1b and c). Firstly, it is necessary to clarify the molecular mechanism of Xinluzhong 54 in response to salt stress. Here, an integrated analysis of transcriptome and proteome was performed to analyze the regulation mechanism of cotton in response to salt stress and identified the key regulatory genes which may be used to improve the salt tolerance in cotton.

Plants initiate defense response against salt stress by regulating the expression of stress responsive genes. It has been widely reported that phytohormones were involved in regulating plant response to salt stress (Ryu and Cho 2015). Transcriptome analysis revealed that cotton response to salt stress induced a variety of hormone signaling pathway genes, reflecting the complexity of hormone involved in the regulation of cotton response to salt stress (Additional file 1: Figure S6). Some evidences suggested that ABA and JA hormone can enhance salt tolerance in many plants (Goossens et al. 2016; Vishwakarma et al. 2017). ABA signaling pathway was activated by salt stress, it was consistent with the increase of ABA content in cotton leaves (Li et al. 2019). In addition, some reports have pointed out that JA also plays an important role in regulating salt tolerance in other plants, such as rice and wheat (Kang et al. 2005; Qiu et al. 2014). A previous study showed that genes involved in JA synthesis and signal transduction were induced by salt stress, such as LOX, AOS, JAZ and MYC. Appropriate concentrations (0.01 μmol·L^− 1) of coronatine, which is structurally and functionally similar to JA, could improve the salt tolerance of cotton (Xie et al. 2015). Overexpression of the JA signal repressor gene GhJAZ2 in cotton significantly reduces the salt tolerance (Sun et al. 2017). The cross regulation between ABA and JA plays an important role in mediating salt tolerance. PnJAZ1, a moss jasmonate ZIM-domain gene is induced by ABA, enhances salt tolerance in plants through inhibiting ABA signaling under salt stress (Liu et al. 2019). The cytokinin regulatory pathway is different from ABA and JA, it has been shown that genes in cytokinin signaling pathways are down-regulated, which are consistent with a previous report (Liu et al. 2012). We speculated that the inhibition of cytokinin signaling may be one of the reasons for repressing cotton growth under salt stress.

Some genes involved in signal transduction were differentially activated or repressed under salt stress, such as transcription factors, CBL interaction protein kinases, MAPK signaling cascade pathway related genes (Fig. 4c and d, Additional file 1: Figure S7). Many transcription factor family genes showed rapid response to salt stress. It was remarkable that most of the differentially expressed NAC genes were up-regulated when cotton encountered salt stress, it is consistent with the previous studies that NAC genes were widely involved in salt response (Sun et al. 2018). Overexpression of salt stress induced NAC gene GHATAF1 enhances the salt tolerance of cotton (He et al. 2016). Transcriptome analysis showed that CBL interaction protein kinase and MAPK signaling genes were induced by salt stress (Fig. 4c, d). In Arabidopsis, CBL2/3 either regulated osmotic or salt stress responses by interacting with CIPK21 and targeting the kinase to the vacuole thus regulating the homeostasis of ions and water on the vacuole (Pandey et al. 2015), or regulated potassium and magnesium ion homeostasis by interacted with CIPK3, CIPK9 and CIPK23 (Liu et al. 2013; Tang et al. 2015). Here, two CBL genes (GhCBL2, GhCBL3), which were homologous to Arabidopsis CBL2 and CBL3, were up-regulated by salt stress and were predicted to interact with several salt-induced CIPK proteins, such as CIPK3, CIPK9 and CIPK23 (Figure S8a). In Arabidopsis, the interaction between MKK2 and MPK4 enhanced salt stress tolerance through MEKK1-MKK2-MPK6/MPK4 pathway (Teige et al. 2004). MAPK cascade MEKK1-MKK2-MPK6/MPK4 could cross-regulate plant response to salt stress with phytohormones, it is involved in SA- and ROS-initiated stress signaling, and plays a vital role in maintaining ROS homeostasis under salt stress (Yang and Guo 2018b). In cotton, the expression patterns of MAPK cascade genes were significantly induced after salt stress treatment at 4 h (Chen et al. 2020). In this study, the homologous genes of Arabidopsis MKK2 and MPK4 genes were up-regulated by salt stress (Fig. 4d). We speculated that some conservative regulatory mechanisms may exist between cotton and Arabidopsis.

Salt stress induced the accumulation of reactive oxygen species and caused lipid peroxidation. In this study, some antioxidant enzyme genes were up-regulated under salt stress (Fig. 5b), which  were consistent with the activation of some antioxidant enzymes in previous studies (Meloni et al. 2003). Salt stress induced the increase of Na⁺ content in cotton leaves (Peng et al. 2014). Here, we identified some transporters and channel proteins encoding genes that showed differently expressed under salt stress, they might play an important role in regulating Na⁺ homeostasis under salt stress (Fig. 5c). Salt stress induced stomatal closure, photosynthetic rate declined, energy metabolism inhibited in plants. Consistent with previous studies, many DEGs and DEPs involved in these biological processes and showed down-regulated expression patterns under salt stress.

Our results showed that only 77, 70 and 47 differentially expressed genes were associated with mRNA and protein levels in T3/T0, T12/T0 and T12/T3, respectively (Fig. 2a). Previous studies have shown that there was a low correlation between transcriptome and proteome results (Chen et al. 2017; Wang et al. 2015). For example, transcriptional factors, CBL-CIPK and MAPK cascade related genes which were identified in transcriptome, were not identified in proteome. The nonlinear and low correlation between these two omics may be due to the following reasons: (i) Post-transcriptional and post-translational modification; (ii) The different rates between transcription and translation, and the different half-life between mRNA and proteins; (iii) The threshold used in the identification of DEGs and DEPs.

However, an integrated analysis of transcriptome and proteome could give more information both on mRNA and protein levels. In the present study, the DEGs which were simultaneously identified at the transcriptional and protein levels were selected for identifying candidate genes. Based on this research strategy, 9 candidate genes that might play an important role in cotton response to salt stress were identified (Fig. 6). AtLEA14 (AT1G01470) is the homologous gene of GhLEA14 (Gh_A11G0835) in Arabidopsis, and the amino acid sequence homology between GhLEA14 and AtLEA14 is 70.86%. Overexpression of AtLEA14 could significantly enhance the salt stress tolerance of Arabidopsis (Jia et al. 2014). We speculate that GhLEA14 is very likely to regulate salt tolerance in cotton. These candidate genes can provide important resources for genetic improvement of salt tolerance in cotton.

Fig. 6

Schematic presentation of a possible regulatory model in cotton (Gossypium hirsutum L. acc. Xinluzhong 54) leaves under salt stress. Cells of cotton sense salt stress by RLKs or other unknown sensors, then transfer the external signal to the cytoplasm, and change the content of Ca²⁺, ROS and other hormones (i.e, ABA, JA). Inside of the cell, Ca²⁺ signal cascade (i.e, GhCBL2,3; GhCIPK8,9,10,25), ROS, MAPK cascade (i.e, GhMAP3K5,13,14,17-GhMKK2,4,9-GhMPK3,4,7) and other hormone signal cascades pathway are activated, which could alter the global transcriptional profiles in cotton (the expression of stress related transcription factors, such as MYB, ERF, NAC and bZIP are initiated). The expression of genes encoding proteins that related to salt stress response, including NHX, osmoregulation related and ROS scavenging proteins (i.e, GhLEA14; GhP5CS; GhPOD; GhSOD; GhCAT; GhFIB; GhADH1). Ultimately, pathways to maintain the osmotic homeostasis, ionic homeostasis and ROS homeostasis are activated to facilitate cotton to resist and adapt to salt stress

Conclusions

According to the integrated analysis of transcriptome and proteome, a complex salt stress response network was outlined in cotton. This study can contribute to understand the mechanism of cotton response to salt stress and facilitate the breeding of salt-resistant cotton varieties by molecular breeding.