The Phosphoproteomic Response of Okra (Abelmoschus esculentus L.) Seedlings to Salt Stress

Soil salinization is a major environmental stresses that seriously threatens land use efficiency and crop yields worldwide. Although the overall response of plants to NaCl has been well studied, the contribution of protein phosphorylation to the detoxification and tolerance of NaCl in okra (Abelmoschus esculentus L.) seedlings is unclear. The molecular bases of okra seedlings’ responses to 300 mM NaCl stress are discussed in this study. Using a combination of affinity enrichment, tandem mass tag (TMT) labeling and high-performance liquid chromatography–tandem mass spectrometry analysis, a large-scale phosphoproteome analysis was performed in okra. A total of 4341 phosphorylation sites were identified on 2550 proteins, of which 3453 sites of 2268 proteins provided quantitative information. We found that 91 sites were upregulated and 307 sites were downregulated in the NaCl/control comparison group. Subsequently, we performed a systematic bioinformatics analysis including gene ontology annotation, domain annotation, subcellular localization, and Kyoto Encyclopedia of Genes and Genomes pathway annotation. The latter revealed that the differentially expressed proteins were most strongly associated with ‘photosynthesis antenna proteins’ and ‘RNA degradation’. These differentially expressed proteins probably play important roles in salt stress responses in okra. The results should help to increase our understanding of the molecular mechanisms of plant post-translational modifications in response to salt stress.


Introduction
Okra (Abelmoschus esculentus L.), an annual herb of Malvaceae family, is native to Africa and India [1,2]. It is not only a nutrient-rich vegetable that is used in traditional Chinese medicines, but also has a high culinary value. As a very important crop and vegetable, it is cultivated in many temperate and subtropical parts of the world [3]. Owing its high oil production rate and great ecological adaptability, okra is a potential bioenergy crop [4]. Most reports focus on its biological characteristics and cultivation techniques [5][6][7], chemical composition and medicinal value [8][9][10], and tissue culture [11,12]. However, there have been few studies on the salt tolerance of okra.
Salt stress is an important environmental condition that limits plant growth and decreases crop productivity. Irrigation water containing trace amounts of sodium chloride (NaCl) can increase soil

Quantitative Phosphoproteomic Data Analysis
Through affinity enrichment and then using the LC-MS/MS approach, the phosphoproteomic changes in okra seedlings treated with salt stress were analyzed. The MS data validation is shown in Figure 1. The mass error of all the identified peptides was confirmed. The distribution of the mass error was near zero and most were less than 0.02 Da, indicating that the mass accuracy of the MS data fitted the requirements ( Figure 1A). The lengths of most peptides were between eight and 20 amino acids (aa), which agreed with the tryptic peptides' properties ( Figure 1B).
The pair-wise Pearson's correlation coefficients among the six samples showed sufficient reproducibility ( Figure S1). The detailed information of identified peptides, including peptide sequences, matching scores, precursor charges, modifications and delta masses, can be found in Table S1. sequences, matching scores, precursor charges, modifications and delta masses, can be found in Table S1.

Annotation and Classification of all Identified Phosphorylated Proteins in Okra
A total of 4341 phosphorylation sites were identified on 2550 proteins, of which 3453 sites of 2268 proteins provided quantitative information (Table S2). Approximately 63.05% of the peptides were modified at a single site, 20.42% were altered at two sites, and 8.54% at three sites ( Figure 1C). To ensure that the results were highly reliable, we filtered the identification data using a criterion of a localization probability >0.75 and ultimately identified 3072 phosphorylation sites on 2175 proteins, of which 2241 sites of 2027 proteins provided quantitative information. The filtered data were used for subsequent bioinformatics analyses. Among the phosphorylation sites, 2607 (84.86%) involved a serine residue, 405 (13.18%) a threonine, and 60 (1.95%) a tyrosine ( Figure 1D).
To understand their functions, all identified phosphorylated proteins were annotated by GO terms based on the three major categories, namely biological process, cellular component and molecular function (Figure 2A). For the biological process category, 'cellular process' (427 proteins), 'metabolic process' (425 proteins), and 'single-organism process' (173 proteins) were the top dominant terms; for the molecular function category, the dominant terms were 'binding' (815 proteins), 'catalytic activity' (458 proteins), and 'transporter activity' (60 proteins); and for the cellular component category, most proteins were related to 'cell' (195 proteins), 'membrane' (149 proteins), and 'organelle' (129 proteins). Most identified proteins were grouped into 15 subcellular component categories predicted by Wolfpsort software, including nucleus (1007 proteins), chloroplast (489 proteins), and cytoplasm (358 proteins) ( Figure 2B). The detailed information, including protein IDs, GO, KEGG, and domain categories, 238 subcellular localizations and functional enrichments of all identified proteins, are listed in Table S3.

Annotation and Classification of all Identified Phosphorylated Proteins in Okra
A total of 4341 phosphorylation sites were identified on 2550 proteins, of which 3453 sites of 2268 proteins provided quantitative information (Table S2). Approximately 63.05% of the peptides were modified at a single site, 20.42% were altered at two sites, and 8.54% at three sites ( Figure 1C). To ensure that the results were highly reliable, we filtered the identification data using a criterion of a localization probability >0.75 and ultimately identified 3072 phosphorylation sites on 2175 proteins, of which 2241 sites of 2027 proteins provided quantitative information. The filtered data were used for subsequent bioinformatics analyses. Among the phosphorylation sites, 2607 (84.86%) involved a serine residue, 405 (13.18%) a threonine, and 60 (1.95%) a tyrosine ( Figure 1D).
To understand their functions, all identified phosphorylated proteins were annotated by GO terms based on the three major categories, namely biological process, cellular component and molecular function (Figure 2A). For the biological process category, 'cellular process' (427 proteins), 'metabolic process' (425 proteins), and 'single-organism process' (173 proteins) were the top dominant terms; for the molecular function category, the dominant terms were 'binding' (815 proteins), 'catalytic activity' (458 proteins), and 'transporter activity' (60 proteins); and for the cellular component category, most proteins were related to 'cell' (195 proteins), 'membrane' (149 proteins), and 'organelle' (129 proteins). Most identified proteins were grouped into 15 subcellular component categories predicted by Wolfpsort software, including nucleus (1007 proteins), chloroplast (489 proteins), and cytoplasm (358 proteins) ( Figure 2B). The detailed information, including protein IDs, GO, KEGG, and domain categories, 238 subcellular localizations and functional enrichments of all identified proteins, are listed in Table S3.

Peptide Motifs Associated with Phosphorylation
A total of 2514 sequences containing the 13 residues, with six upstream and six downstream residues surrounding each of the phosphorylation sites, were obtained ( Figure S2). Of these, 2255 (89.70%) were centered on a serine residue, 237 (9.43%) on a threonine, and 22 (0.87%) on a tyrosine. The serine phosphorylation category included 29 overpresented motifs: the most frequent motifs were 'sP' (355 occurrences) and 'PxsP' (182 occurrences), followed by 'Rxxs' (151 occurrences) and 'RSxs'(116 occurrences). 'tP' was the most frequent motif in the threonine phosphorylation category.

Differentially Phosphorylated Proteins in Response to NaCl Treatment
To compare the differentially phosphorylated proteins (DPPs) between different samples, expression profiles of the proteins are shown in a heatmap ( Figure 3A). To reflect the changing trends among the six samples, all of the DPPs were categorized into four clusters with MeV software using the K-means method. The proteins in cluster I showed high accumulations in the NaCl stress-treated samples, while the rest of the proteins were downregulated in the NaCl stress-treated samples ( Figure 3B). Among the DPPs, 72 proteins were upregulated and 270 proteins were downregulated at 48 h after NaCl treatment compared with the control ( Figure 3C). All DPPs were classified into 11 subcellular components according to their subcellular localizations ( Figure 2B), including 155 nuclear-localized DPPs, 86 chloroplast-localized DPPs, and 47 cytoplasm-localized DPPs. The detail information of DPPs and differential phosphosites are listed in Table S4

Peptide Motifs Associated with Phosphorylation
A total of 2514 sequences containing the 13 residues, with six upstream and six downstream residues surrounding each of the phosphorylation sites, were obtained ( Figure S2). Of these, 2255 (89.70%) were centered on a serine residue, 237 (9.43%) on a threonine, and 22 (0.87%) on a tyrosine. The serine phosphorylation category included 29 overpresented motifs: the most frequent motifs were 'sP' (355 occurrences) and 'PxsP' (182 occurrences), followed by 'Rxxs' (151 occurrences) and 'RSxs'(116 occurrences). 'tP' was the most frequent motif in the threonine phosphorylation category.

Differentially Phosphorylated Proteins in Response to NaCl Treatment
To compare the differentially phosphorylated proteins (DPPs) between different samples, expression profiles of the proteins are shown in a heatmap ( Figure 3A). To reflect the changing trends among the six samples, all of the DPPs were categorized into four clusters with MeV software using the K-means method. The proteins in cluster I showed high accumulations in the NaCl stress-treated samples, while the rest of the proteins were downregulated in the NaCl stress-treated samples ( Figure 3B). Among the DPPs, 72 proteins were upregulated and 270 proteins were downregulated at 48 h after NaCl treatment compared with the control ( Figure 3C). All DPPs were classified into 11 subcellular components according to their subcellular localizations ( Figure 2B), including 155 nuclear-localized DPPs, 86 chloroplast-localized DPPs, and 47 cytoplasm-localized DPPs. The detail information of DPPs and differential phosphosites are listed in Table S4. To predict the biological functions, all DPPs were assigned to GO terms. For upregulated proteins, the largest numbers of DPPs were found in 'metabolic process' (eight proteins) in the biological process GO term, in cellular component, the most frequent was 'membrane' (six proteins), while, in molecular function, 'binding' (19 proteins) was the most highly represented group ( Figure  4A). For downregulated proteins, the most enriched biological process GO term was 'cellular process'; within "cellular component", the most enriched was related to 'cell' (36 proteins); and the most enriched term within molecular function was 'binding' (91 proteins) ( Figure 4B). In total, 41 DPPs were involved in membrane and transport proteins according to the GO terms and their subcellular localizations ( Table 1). Eight of them were induced by the NaCl treatment. The GO enrichment analysis revealed that the molecular function GO term of upregulated DPPs most significantly enriched was 'oxidoreductase activity, acting on NAD(P)H' and that the biological process' GO term most significantly enriched was 'metal ion transport' ( Figure S3). For the downregulated DPPs, the molecular function GO terms most significantly enriched was 'transferase activity, transferring glycosyl groups' while the biological process GO term most significantly enriched was 'DNA metabolic process' ( Figure S3).  To predict the biological functions, all DPPs were assigned to GO terms. For upregulated proteins, the largest numbers of DPPs were found in 'metabolic process' (eight proteins) in the biological process GO term, in cellular component, the most frequent was 'membrane' (six proteins), while, in molecular function, 'binding' (19 proteins) was the most highly represented group ( Figure 4A). For downregulated proteins, the most enriched biological process GO term was 'cellular process'; within "cellular component", the most enriched was related to 'cell' (36 proteins); and the most enriched term within molecular function was 'binding' (91 proteins) ( Figure 4B). In total, 41 DPPs were involved in membrane and transport proteins according to the GO terms and their subcellular localizations ( Table 1). Eight of them were induced by the NaCl treatment. The GO enrichment analysis revealed that the molecular function GO term of upregulated DPPs most significantly enriched was 'oxidoreductase activity, acting on NAD(P)H' and that the biological process' GO term most significantly enriched was 'metal ion transport' ( Figure S3). For the downregulated DPPs, the molecular function GO terms most significantly enriched was 'transferase activity, transferring glycosyl groups' while the biological process GO term most significantly enriched was 'DNA metabolic process' ( Figure S3).

Enrichment Analysis of the DPPs under the NaCl Treatment
The KEGG enrichment analysis revealed that the upregulated DPPs were most strongly associated with 'Photosynthesis-antenna proteins' and 'Starch and sucrose metabolism' pathways. The downregulated DPPs were most strongly associated with 'Homologous recombination' and 'RNA degradation' pathways ( Figure 5A). We further summarized the DPPs involved in photosynthesis and RNA degradation ( Figure 6). In total, six proteins related to photosynthesis were identified, including light-harvesting complex II (LHCII) chlorophyll a/b-binding protein 1/2/4, photosystem II oxygen-evolving enhancer protein 2, photosystem II protein H (PsbH), and PSI

Enrichment Analysis of the DPPs under the NaCl Treatment
The KEGG enrichment analysis revealed that the upregulated DPPs were most strongly associated with 'Photosynthesis-antenna proteins' and 'Starch and sucrose metabolism' pathways. The downregulated DPPs were most strongly associated with 'Homologous recombination' and 'RNA degradation' pathways ( Figure 5A). We further summarized the DPPs involved in photosynthesis and RNA degradation ( Figure 6). In total, six proteins related to photosynthesis were identified, including light-harvesting complex II (LHCII) chlorophyll a/b-binding protein 1/2/4, photosystem II oxygen-evolving enhancer protein 2, photosystem II protein H (PsbH), and PSI reaction center subunit II. In total, 11 phosphorylated sites were identified in these proteins ( Figure 6A). We also measured some physiological parameters owing to differences in some proteins in the photosynthetic pathways. The Fv/Fm image is shown in Figure 6B. The maximum quantum efficiency of photosystem II (PSII) (Fv/Fm) and Pn values decreased significantly in NaCl treatments compared with their respective controls, but the SPAD values did not change significantly. Five proteins related to RNA degradation were identified, including DNA topoisomerase 2-associated protein PAT1, mRNA-decapping enzyme subunit 2 (DCP), enhancer of mRNA-decapping protein 4, CCR4-NOT transcription complex subunit 2, and CCR4-NOT complex subunit 16 ( Figure 6C). reaction center subunit II. In total, 11 phosphorylated sites were identified in these proteins ( Figure  6A). We also measured some physiological parameters owing to differences in some proteins in the photosynthetic pathways. The Fv/Fm image is shown in Figure 6B. The maximum quantum efficiency of photosystem II (PSII) (Fv/Fm) and Pn values decreased significantly in NaCl treatments compared with their respective controls, but the SPAD values did not change significantly. Five proteins related to RNA degradation were identified, including DNA topoisomerase 2-associated protein PAT1, mRNA-decapping enzyme subunit 2 (DCP), enhancer of mRNA-decapping protein 4, CCR4-NOT transcription complex subunit 2, and CCR4-NOT complex subunit 16 ( Figure 6C).

Protein-Protein Interactions (PPIs) of DPPs
A PPI network analysis can reveal the relationship between biological functions of different phosphorylated proteins. The high quality PPI maps are shown in Figure 7, and the detailed node and network information are listed in Tables S5 and S6. A total of 39 DPPs (76 modified sites), including five up-and 34 downregulated peptides, was displayed in the PPI network (Figure7). Three enriched interaction clusters were identified from our analysis. For cluster 1, six 'phosphorylation'-related proteins (downregulated), three 'binding'-related proteins (downregulated), a small subunit ribosomal protein S6e (downregulated), a protein-tyrosine phosphatase-like protein, a WD40-repeat-containing domain protein (downregulated), and a PPM-type phosphatase domain (upregulated), have been included. Four 'Ribosome'-related proteins (downregulated) and three 'Spliceosome'-related proteins (upregulated) were identified in cluster 2. For cluster 3, six nucleus-localized proteins, three mitochondria-localized proteins, two chloroplast-localized proteins, and a chloroplast-localized protein cytoplasm have been identified. All DPPs in cluster three were downregulated.

Protein-Protein Interactions (PPIs) of DPPs
A PPI network analysis can reveal the relationship between biological functions of different phosphorylated proteins. The high quality PPI maps are shown in Figure 7, and the detailed node and network information are listed in Tables S5 and S6. A total of 39 DPPs (76 modified sites), including five up-and 34 downregulated peptides, was displayed in the PPI network ( Figure 7). Three enriched interaction clusters were identified from our analysis. For cluster 1, six 'phosphorylation'-related proteins (downregulated), three 'binding'-related proteins (downregulated), a small subunit ribosomal protein S6e (downregulated), a protein-tyrosine phosphatase-like protein, a WD40-repeat-containing domain protein (downregulated), and a PPM-type phosphatase domain (upregulated), have been included. Four 'Ribosome'-related proteins (downregulated) and three 'Spliceosome'-related proteins (upregulated) were identified in cluster 2. For cluster 3, six nucleus-localized proteins, three mitochondria-localized proteins, two chloroplast-localized proteins, and a chloroplast-localized protein cytoplasm have been identified. All DPPs in cluster three were downregulated.

Discussion
Plant responses to salt stress are mediated by complex molecular mechanisms, including signal transduction, and the transcription and translation of salt stress-related genes. Reversible protein phosphorylation is a key regulatory mechanism operating in response to abiotic and biotic stresses in plant [40,41]. In this study, a TMT-based quantitative phosphoproteomic approach was used to study the response of okra shoots to salt stress (Figure 1). Many phosphorylation sites were identified, and some phosphorylated proteins that might be involved in salt stress response were identified (Figure 3). These salt-responsive phosphoproteins may play important roles in salt stress signaling and response in okra shoots.
The first method for large-scale phosphorylation proteomics involved two-dimensional gel electrophoresis, with spots identified by mass spectrometry [28]. However, it is difficult to identify proteins of low abundance, low molecular weight (<15 kDa) or high molecular weight (>150 kDa); superacidic or basic proteins; or hydrophobic proteins by this method [42]. Recently, an MS/MS-based analysis strategy using isobaric tag for relative and absolute quantitation or TMT-labeling was developed for large-scale phosphorylated protein quantification [43,44]. In our study, 4341 phosphorylation sites, 2550 phosphorylated proteins, and 342 DPPs were identified

Discussion
Plant responses to salt stress are mediated by complex molecular mechanisms, including signal transduction, and the transcription and translation of salt stress-related genes. Reversible protein phosphorylation is a key regulatory mechanism operating in response to abiotic and biotic stresses in plant [40,41]. In this study, a TMT-based quantitative phosphoproteomic approach was used to study the response of okra shoots to salt stress (Figure 1). Many phosphorylation sites were identified, and some phosphorylated proteins that might be involved in salt stress response were identified (Figure 3). These salt-responsive phosphoproteins may play important roles in salt stress signaling and response in okra shoots.
The first method for large-scale phosphorylation proteomics involved two-dimensional gel electrophoresis, with spots identified by mass spectrometry [28]. However, it is difficult to identify proteins of low abundance, low molecular weight (<15 kDa) or high molecular weight (>150 kDa); superacidic or basic proteins; or hydrophobic proteins by this method [42]. Recently, an MS/MS-based analysis strategy using isobaric tag for relative and absolute quantitation or TMT-labeling was developed for large-scale phosphorylated protein quantification [43,44]. In our study, 4341 phosphorylation sites, 2550 phosphorylated proteins, and 342 DPPs were identified (Table S2). The large number of identified phosphorylated proteins gave us the opportunity to conduct a more in-depth and comprehensive analysis of salt-responsive proteins. 'sP' and 'Rxxs' are the two most frequently occurring motifs for phosphoserine, and this was confirmed in okra. The 'sP' motif, which provides a target for cAMP-and cGMP-dependent protein kinase C, CDK protein kinase 2, Ca-dependent protein kinases, mitogen-activated protein kinase, receptor-like kinase, sucrose nonfermenting 1-related protein kinase 2 and STE20-like kinase, has also been identified as being overrepresented in other plants [45] (Figure S2). Meanwhile, the 'Rxxs' motif has been reported to be associated with mitogen-activated protein kinase kinase, protein kinase A, and calmodulin-dependent protein kinase [35,45].
In the face of various biotic and abiotic stresses, plants gradually regulate water and ion transport mechanisms [46]. Maintaining cell osmotic potential under salt stress is a major challenge for plant growth and development [47]. In this study, several membrane proteins and ion transporters were shown to be differentially phosphorylated under salt stress, suggesting their phosphorylation of these proteins may regulate osmotic balance in okra under salt stress (Table 1). Under salt stress, several transporters of sugars, amino acids or ions were differentially phosphorylated in okra, and some might be involved in osmotic regulation. Three phosphorylation sites of two multidrug/pheromone exporter proteins in the ATP-binding cassette (ABC) transporter subfamily B decreased, and a white-brown complex ABC transporter protein was induced after the salt treatment in this study. ABC transporters are a class of transmembrane proteins that are ubiquitous in prokaryotes and eukaryotes [48]. Most ABC transporters bind to adenosine triphosphate (ATP) and release energy by hydrolyzing ATP to regulate the transmembrane transport of substances [48]. Plant ABC transporters play key roles in plant secondary metabolite transport and accumulation, phytohormone transport, lipid metabolism, exogenous toxin detoxification, and plant disease resistance [49]. Some transporters and iron pumps, such as an amino acid transporter (Unigene80449_All), a sugar phosphate transporter (Unigene66780_All), an oligopeptide transporter (Unigene78857_All), mechanosensitive ion channel protein 8-like (CL18.Contig22_All), phosphate transporter PHO1 homolog 3(CL27908.Contig8_All), a sugar transporter (CL2495.Contig20_All), UDP-glucuronic acid/UDP-N-acetylgalactosamine transporter (CL7033.Contig3_All), and vacuolar cation/proton exchanger 1a (CL27084.Contig3_All) were repressed by salt stress. Changes in phosphorylation status of these transporters after salt stress treatment suggest that they may play a role in regulating ion and small molecule concentrations to balance cell osmotic potential and reduce cytoplasmic water loss. Four phosphorylation sites of three cellulose synthase proteins decreased after salt treatment in this study. The cell wall plays an important role in maintaining cell morphology, and the cellulose synthase family proteins is involved in cell wall formation was regulated at phosphorylation level [50,51].
We detected changes in some phosphorylated proteins involved in photosynthesis after the salt treatment ( Figure 6A). Similar results were also reported found in A. thaliana [24]. Phosphorylation at sites on three light-harvesting II complex (LHCII) proteins-LHB1B2, LHCb4.2, and LHCb1.2-were induced following exposure to salt or H 2 O 2 stresses, as were sites on PsbH and PsaP [24]. The phosphorylation of LHCII proteins adapts to environmental changes [52]. The phosphorylation of LHCII participates in the energy distribution between the two photosystem, and the signal transduction between light reception and the phosphorylation of LHCII is correlated with the redox state of the proton quinone pool. Thus, a reduction of proton quinone leads to the activation of the kinase [53]. Phosphorylation of PsbH proteins regulates the whole membrane network rearrangement of plant thylakoids [24]. Therefore, the phosphorylation of chloroplast proteins affects the redox state [54]. Several DDPs involved in the mRNA degradation pathway ( Figure 6C). Removing the 5 caps of mRNA is an important step in post-transcriptional regulation, because it represents a movement away from active translation and is a prerequisite for the rate-limiting degradation of exogenous ribonucleases [55]. The correlation between transcriptional level and protein abundance was poor after the salt treatment [56]. The transcriptional level was mainly disconnected from the active polymer during dehydration stress [57], indicating that the post-transcriptional regulation mechanism was activated in the early dehydration reaction [58,59]. Our data identified that phosphorylation of DNA topoisomerase 2-associated protein PAT1 increased at the threonine residue, while mRNA-decapping enzyme subunit 2, enhancer of mRNA-decapping protein 4, CCR4-NOT transcription complex subunit 2, and CCR4-NOT complex subunit 16 decreased at serine residues following the NaCl treatment, which suggests that mRNA decapping is altered in response to salt stress.
Ca 2+ is involved in plant signal transduction under environmental stress, and CDPks play a key role in Ca 2+ mediated signal transduction pathway [51]. A calmodulin-like protein(CL9944.Contig4_All) and CDPk protein (Unigene71267_Alll) were decreased at phosphorylation level under NaCl treatments, indicating the Ca 2+ -mediated pathway participates in salt stress signaling in okra shoots ( Figure 7 and Table S5). Protein degradation in plants may be activated under salt stress. Previous proteomic studies have identified a large number of salt stress-related proteins related to protein synthesis, folding and degradation [47,60,61]. In the present study, we observed that five ribosomal protein (40S ribosomal protein S6A, CL2686.Contig12_All; ribosomal protein S15e, CL9991.Contig2_All; ribosomal protein S9, CL19354.Contig5_All; ribosomal protein S6e, Unigene79731_All; ribosomal protein S6-like, Unigene55159_All) decreased at the phosphorylation level under NaCl treatment ( Figure 7 and Table S4). The reduction of these phosphoribosomal proteins supports the view that salt stress usually inhibits protein synthesis [36]. Ribosomal protein S6 (rps6), a component of the 40S ribosomal subunit, whose phosphorylation status can be affected by osmotic stress via TOR pathway [62,63].

Plant Materials and Treatment
Seeds of the okra cultivar 'Wufu' were used in this study. The seeds were disinfected with 1% sodium hypochlorite for 10 min and washed with distilled water three times. Seeds were sown in plastic trays of turfy soil. One week after germination, half-strength Hoagland's nutrient solution was used for irrigation every 3 days. Seedlings were grown in an artificial illumination incubator with a photoperiod of 12-h light/12-h dark, relative humidity of 60%, and light intensity of 300 µmol m −2 s −1 . Three weeks after germination, 20 mL of 300 mmol L −1 NaCl was applied. After 48 h treatment, the chlorophyll fluorescence images were captured using a FluorCam 800 imaging system (PhotonSystems Instruments, Brno, Czechia Republic) as previously described by Feng et al. [64]. After acquiring the chlorophyll fluorescence imaging data, the aboveground parts of the okra plants were harvested to measure the chlorophyll content using a SPAD meter (502DL Plus, SPECTRUM, Illinois, IL, USA). Net photosynthesis (Pn) was measured by using a portable photosynthesis system (LI 6400, LI-COR, Lincoln, NE, USA) with a red/blue LED light source at 1000 µmol m −2 s −1 Photosynthetic Photon Flux Density.

Protein Extraction and Digestion
At 48 h after salt treatment, the aboveground seedling samples were fully ground to a powder in liquid nitrogen. Four-fold volume phenol extraction buffer (containing 10 mM dithiothreitol, 1% protease inhibitor, and 1% phosphatase inhibitor) was added to the samples of each group, and ultrasonic pyrolysis was performed in high intensity ultrasonic processor (Scientz, Ningbo, China). An equal volume of Tris-saturated phenol was added and centrifuged (4 • C, 10 min, 5000 g). The supernatant was removed and precipitated overnight with a 5× volume of 0.1 M ammonium acetate/methanol solution. The precipitated protein was washed successively with methanol and acetone. The protein was redissolved in 8 M urea, and the protein concentration was determined using a bicinchoninic acid assay (BCA) kit (CW0014; CWBIO, Beijing, China), according to the manufacturer's instructions.
For trypsin digestion, dithiothreitol was added to the protein solution, resulting in a 5-mM final concentration, and reduced for 30 min at 56 • C. Then, iodoacetamide was added to an 11-mM final concentration and the samples were incubated at room temperature for 15 min in the dark. Finally, the urea concentration of the samples was diluted to below 2 M by adding 100 mM triethyl ammonium bicarbonate. Trypsin was added at 1:50 mass ratio (trypsin: protein) for the first hydrolysis overnight at 37 • C and at 1:100 (trypsin: protein) mass ratio for a second 4-h digestion.

TMT Labeling and High-Performance Liquid Chromatography (HPLC) Fractionation
After trypsin digestion, peptides were desalted using Strata X C18 solid-phase extraction (Phenomenex, Torrance, CA, USA) and vacuum freeze-dried. The peptides were dissolved in 0.5 M triethyl ammonium bicarbonate and labeled using a TMT kit (ThermoFisher, Shanghai, China) according to the manufacturer's instructions. Briefly, the labeled reagent was dissolved in acetonitrile after thawing and incubated at room temperature for 2 h after mixing with the peptides. Then, the peptide mixture was desalted and freeze-dried under vacuum.
The peptides were fractionated by high pH reversed HPLC with a Thermo Betasil C18 column (5 µm diameter, 10 mm inner diameter, 250 mm length). The operation was as follows: peptides were first separated into 60 fractions using a gradient of 8% to 32% acetonitrile (pH 9.0) over 60 min. Then, the peptides were pooled into 18 fractions and dried by vacuum centrifugation.

Affinity Enrichment
The peptides were dissolved in a concentrated buffer solution (50% acetonitrile/6% trifluoroacetic acid), and the supernatant was transferred to the immobilized metal ion/metal chelate material, which had been washed in advance. The peptide-resin mixture were placed on a rotary table and gently shaken. After this incubation period, the resin was washed three times successively with buffer solutions of 50% acetonitrile/6% trifluoroacetic acid and 30% acetonitrile/0.1% trifluoroacetic acid. The phosphopeptides were then eluted with 10% ammonia solution. The eluent was collected and vacuum freeze-dried. After being drained, the eluent was de-salted using C18 ZipTips according to the manufacturer's instructions. After vacuum freeze-drying, the eluent was used for liquid chromatography-tandem mass spectrometry (LC-MS/MS) analysis.

LC-MS/MS Analysis
The peptides were dissolved in an aqueous solution of 0.1% (v/v) formic acid and separated using an EASY-nLC 1200 ultrahigh-performance liquid phase system (Thermo Fisher Scientific, Waltham, MA, USA). The liquid gradient settings were as follows: 0-38 min, 4%-22% B; 38-52 min, 22%-32% B; 52-56 min, 32%-80% B. Solvent B was an aqueous solution containing 0.1% formic acid and 90% acetonitrile. The flow rate was maintained at 450 nL/min. The peptides were ionized using an NSI ion source and then analyzed by Q Exactive TM HF-X mass spectrometry (ThermoFinnigan, Somerset, NJ, USA). The ion source voltage was set at 2.0 kV, and the peptide parent ions and their secondary fragments were detected and analyzed using a high-resolution Orbitrap. The scanning range of the first-order mass spectrometry was set to 350-1600 m/z, and the scanning resolution was set to 60,000. The fixed starting point of the second-order mass spectrum scanning range was 100 m/z, and the second-order scanning resolution was set at 30,000. In the data acquisition mode, the data-dependent scanning program was used, in which, the parent ions of the first 20 peptides with the highest signal intensity were selected to enter the higher-energy C-trap dissociation collision pool successively after the first scan. The fragmentation energy was set at 28%, and the second-order mass spectrometry analysis was performed successively. To improve the effective utilization rate of mass spectrometry, the automatic gain control was set to 1E5, the signal threshold was set to 5E3, the maximum injection time was set to 100 ms, and the dynamic exclusion time of the tandem mass spectrometry scanning was set to 15 s to avoid the repeated scanning of parent ions.

Database Search
The resulting MS/MS data were queried against our previously published okra transcriptome data (NCBI Sequence Read Archive database accession: SRP130180) using the Maxquant search engine (v.1.5.2.8) concatenated with the reverse decoy database [65]. The enzyme digestion mode was set as Trypsin/P. The number of missing cut points was set as 2. The minimum length of the peptide segment was set to seven amino acid (aa) residues. The maximum modification number of modifications for a peptide segment was set as five. First search and Main search were set as 20 PPM and 5 PPM, respectively, and the mass error tolerance of secondary debris ions was 0.02 Da. Cysteine alkylation was set as the fixed modification, and the alterable modification was methionine oxidation. The quantitative method was set as tmt-6plex, and the FDR for protein identification and PSM identification were set at 1%.

Annotation Methods and Functional Enrichment
The gene ontology (GO) annotation of the proteome was derived from the UniProt-GOA database (http://www.ebi.ac.uk/GOA/). First, the identified protein IDs were converted to UniProt IDs and then mapped to GO IDs. The InterProScan software was used to annotate the protein's GO function based on the protein sequence alignment method if the proteins were not annotated by the UniProt-GOA database. For each category, a two-tailed Fisher's exact test was employed to test the enrichment of the differentially expressed protein (DEP) against all identified proteins. A GO annotation with a corrected p-value < 0.05 was considered significant. The Kyoto Encyclopedia of Genes and Genomes (KEGG) database was used to annotate protein pathways. The InterPro domain database (http://www.ebi.ac.uk/interpro/) was used to analyze identified the functional descriptions of protein domains. A pathway with a corrected p-value < 0.05 was considered significant. Wolfpsort (http://www.genscript.com/psort/wolf_psort.html), a subcellular localization predication software, PSORT/PSORT II version, was used to predict subcellular localization. The InterPro database was searched, and a two-tailed Fisher's exact test was employed to test the enrichment of the DEP against all identified proteins. Protein domains with a p-value < 0.05 were considered significant. The software motif-x was used to analyze the models of sequences that contained the amino acids in specific positions of modify-13-mers (six amino acids upstream and downstream of the site) in all protein sequences. Additionally, all the database protein sequences were used as background database parameters, and the other parameters were set at default.

Protein-Protein Interaction Network
All differentially phosphorylated proteins that were classified as identifiers were queried against the STRING database version 10.5 (http://string-db.org/cgi/input.pl) for protein-protein interactions. The differential protein interactions were extracted according to the confidence score of >0.7 (high confidence). Then, the R package "network workd3" tool (https://cran.r-project.org/web/packages/ networkD3/) was used to visualize the differential protein interaction networks.

Conclusions
In summary, we presented a comprehensive large-scale phosphoproteome of okra using a combination of affinity enrichment and TMT labeling high-resolution LC-MS/MS analysis. A total of 4341 phosphorylation sites were identified on 2550 proteins. We also investigated the DPPs that changed after the exposure of okra seedlings to salt or water for 48 h. In total, 342 DPPs (399 phosphorylated sites) were identified, 72 (91 sites) of which were upregulated and 270 (307 sites) were downregulated under salt stress conditions. A number of DPPs were involved in photosynthesis and mRNA degradation pathway. Our data will provide a basic resource that will assist understanding of the molecular mechanisms involved in the responses of okra plants to salt stress.