Dynamic Phosphoproteome Analysis of Seedling Leaves in Brachypodium distachyon L. Reveals Central Phosphorylated Proteins Involved in the Drought Stress Response

Drought stress is a major abiotic stress affecting plant growth and development. In this study, we performed the first dynamic phosphoproteome analysis of Brachypodium distachyon L. seedling leaves under drought stress for different times. A total of 4924 phosphopeptides, contained 6362 phosphosites belonging to 2748 phosphoproteins. Rigorous standards were imposed to screen 484 phosphorylation sites, representing 442 unique phosphoproteins. Comparative analyses revealed significant changes in phosphorylation levels at 0, 6, and 24 h under drought stress. The most phosphorylated proteins and the highest phosphorylation level occurred at 6 h. Venn analysis showed that the up-regulated phosphopeptides at 6 h were almost two-fold those at 24 h. Motif-X analysis identified the six motifs: [sP], [Rxxs], [LxRxxs], [sxD], [sF], and [TP], among which [LxRxxs] was also previously identified in B. distachyon. Results from molecular function and protein-protein interaction analyses suggested that phosphoproteins mainly participate in signal transduction, gene expression, drought response and defense, photosynthesis and energy metabolism, and material transmembrane transport. These phosphoproteins, which showed significant changes in phosphorylation levels, play important roles in signal transduction and material transmembrane transport in response to drought conditions. Our results provide new insights into the molecular mechanism of this plant’s abiotic stress response through phosphorylation modification.


Results
Phenotypic, physiological, and ultrastructural changes in seedling leaves under drought conditions. The dynamic changes in phenotypes and some important physiological indices of Brachypodium distachyon 21 (Bd21) seedling leaves under different drought conditions are shown in Figs 1 and S1. With the intensification of drought stress, leaves and roots displayed more obvious changes. At 24 h, leaves had a slight curl, while roots became thin compared to the control group. There were no significant differences at 6 h (Fig. 1A). Drought stress resulted in significant changes in the physiological index of seedling leaves (Figs S1 and 1B). The results from five different periods (0, 6, 12, 24, and 48 h) showed that relative water content (RWC) and total chlorophyll content for the 20% PEG 6000 treatment displayed a continuous declining trend. In contrast, the levels of malondialdehyde (MDA) and proline (Pro) content continuously increased in response to drought stress. The 24 h treatment was a critical turning point, while 6 h was a transition period. Compared to the control group (0 h), MDA exhibited a significant (p < 0.05) and highly significant (p < 0.01) increase at 6 h (1.299-fold) and 24 h (1.939-fold), respectively. Proline content was significantly increased by 1.646-fold at 24 h (p < 0.05).
To reduce experimental error, unified flag leaves of B. distachyon were selected for transmission electron microscope (TEM) observation. The ultrastructural changes between the non-stress and drought stress treatments of Bd21 seedling leaves are shown in Fig. 1C. The flag leaves in the control group had normal whole mesophyll cells and clear subcellular structures. Treatment with 20% PEG resulted in different degrees of damage to parts of the mesophyll cells and some organelles, including clear plasmolysis of the cell wall, misshapen chloroplasts and nucleolus, particularly at 24 h (Fig. 1C). With disrupted thylakoids, the number of starch granules increased, mitochondria became swollen, and the color faded. These results demonstrate that drought stress has significant effects on leaf ultrastructure.
To increase reliability and accuracy, only phosphorylated sites belonging to the highest phosphorylation site localization probability (p ≥ 0.9, class I) were used for further analysis. Phosphopeptides with receivable Student's t test (p < 0.05) consequence, significant intensity changes (≥ 2-fold), high phosphorylation site localization probability (≥ 0.9), and phosphorylation site score differences exceeding five were taken into consideration to identify significant changes in phosphorylation level. In the strict limitations screening progress, a total of 484 phosphorylation sites had significant changes in the phosphorylation level (SCPL), and divided into 0-6 h and 6-24 h. Comparisons of phosphorylation sites, phosphopeptides, and phosphorylated proteins are provided in Fig. 2a~c, Tables S2a and S2b, respectively. In two groups, there were 76 phosphorylation sites in common, representing 82 phosphorylated proteins.
Scientific RepoRts | 6:35280 | DOI: 10.1038/srep35280 A Venn diagram was used to analyze the significantly up-regulated and down-regulated phosphopeptides in both groups (Fig. 2b). The fold change in up-regulated phosphopeptides at 6 h was almost two-fold greater than at 24 h. Six and two phosphopeptides were up-regulated and down-regulated, respectively (Fig. 2d). Proteins involved in stress defense included a stress protein DDR48-like (gi|357163456) and a late embryogenesis abundant protein B19.3 (LEA, gi|721639026).
Expression pattern analysis of phosphoproteins with SCPL under drought stress. Principal component analysis (PCA) is a multivariate method of dimensionality that can be used to identify relationships among variables 25 . The PCA results from 442 SCPL phosphoproteins are shown in Fig. 3. Three groups of biological replication (NP1, NP2, and NP3) showed a high level of reproducibility (Fig. 3a). The treatment groups were significantly different from the control group. The treatment group at 6 h had the greatest differences from the other treatments. Some spots focused on the origin of coordinates and a few had scattered distribution in the periphery (Fig. 3b). Most of the scattered spots belonged to pattern I and pattern II in Fig. S5. Meanwhile, these scattered spots were common spots of 76 phosphorylation sites. These results indicate that Bd21 underwent significant changes in protein phosphorylation levels under drought stress. To visualize the dynamic changes in the SCPL phosphopeptides, the Euclidean distance similarity metric was used to define the similarity, and hierarchical clusters were assembled using the complete linkage clustering method. Five hierarchical clusters (patterns I-V), corresponding to the different expression patterns, were defined in each period ( Fig. S5; Table S3). Pattern I, containing 179 phosphorylation sites, was significantly up-regulated at 6 h, while pattern II, with 106 phosphorylation sites, showed up-regulated expression at 24 h. Furthermore, we classified the function of these two parts phosphorylated proteins in Figs S7a,b. These results further validate the previous PCA analysis, where the 6 h treatment had significantly more SCPL proteins compared to the other treatments.
Phosphopeptide screening and placement in functional categories. All of the identified SCLP proteins were used to predict the functional categories with GO annotation. They were involved in different biological processes, molecular functions, and cellular components ( Fig. 4; Table S4). Regarding biological processes, primary metabolic processes accounted for 20% in Fig. 4a. The proportion of gene expression and single-organism cellular process was less than primary metabolic process, respectively. Nucleotide and nucleic  acid binding occupied the main part of the molecular functions (Fig. 4b), indicating that nucleoproteins respond to drought stress positively and phosphorylation is one of the coping mechanisms. The cellular components of more than 30% of the phosphoproteins were localized to the membrane (Fig. 4c), indicating that they play a key role in metabolism and some phosphoproteins of water and ion transport were identified in Table 1. Furthermore, three aquaporins (AQPs) were analyzed by 3D structural analyses and sequence alignment in Fig. S6. These AQPs were aquaporin PIP1-5 (IPI1-5, gi|357137703), probable aquaporin IPIP2-7 (IPIP2-7, gi|357159722), and aquaporin NIP2-2 (NIP2-2, gi|357124731). In 3D structural analyses, only PIP1-5 phosphorylation site belonged to the α -helix. Others were on the random coil, which is in agreement with the secondary structure surrounding the phosphorylation site prediction results (Fig. 5b). Most of the SCPL phosphoproteins (nearly 90%) were random coils or turns, 7.6% were α -helixes, and the remaining proteins (1.6%) were β -pleated sheets ( Fig. 5b; Table S5b). Some SCPL proteins identified by Lv et al. 26 under salt stress were also identified under drought stress in this study ( Table 2). Their functional classification is provided in Fig. S7c. A total of 51 phosphorylated proteins were identified; 31% of them are related to transcription and translation, and approximately 18% participate in water and ion transport related to osmotic regulation. However, 29% of the proteins are associated with other processes and require further study (Fig. S7c).
Phosphorylation motif analysis of significantly changed phosphopeptides. Phosphorylation motifs of the phosphoproteins linked to kinase were ascertained by Motif-X and WebLogo online software tools ( Fig. 5a; Table S5a). There were six motifs were identified, including five phosphorylation residues were enriched by Ser, ( 28 . In addition, motif [sxD] was an acidic motif recognized by casein kinase-II (CKII) 29 . Only motif [sF] was easily recognized by protein kinase C 30 .
Conservation analysis of phosphoproteins responsive to drought stress. Sequences of the SCPL phosphoproteins were used as queries to BLAST search the phosphoprotein databases constructed based on datasets in the Plant Protein Phosphorylation DataBase (P3DB) 31 , Medicago-Omics Repository (MORE) 32 , and PhosPhAt 4.0 33 . Oryza. sativa, A. thaliana, Medicago truncatula, and others were compared to B. distachyon to determine the extent of conservation of phosphoproteins among different plant species. The thresholds were set at a score ≥ 80, E-value < 1E− 10, and identity ≥ 30% ( Fig. 5c; Table S5c). In total, among the 442 phosphoproteins identified in B. distachyon, 104 (23.48%), 86 (19.41%), and 49 (11.06%) had phosphorylated orthologs in three, two, and one species, respectively. Only 164 phosphoproteins had phosphorylated orthologs in all four species. A total of 40 (9.02%) phosphoproteins had no phosphorylated orthologs in any of these species. Protein-protein interaction (PPI) analysis. To explore the relationships in the PPI network, 442 SCPL proteins identified in the current study were analyzed by STRING. A total of 307 EuKaryotic Orthologous Groups (KOGs) were used to construct the PPI network. To improve the reliability of PPI analysis, the confidence score was set at the highest level (≥ 0.900). Finally, a complex PPI network was completed through Cytoscape ( Fig. S8; Table S6).
Numerous SCPL phosphoproteins were found to participate in transcription and translation of gene expression (Fig. S8). Some of the phosphoproteins belonged to nucleotide or nucleic acid binding phosphorylation proteins, in agreement with the functional categories ( Fig. 4b). Moreover, some of the important protein kinases and phosphatases, including abscisic acid-inducible protein kinase (PKAB1, KOG0583), myotonin protein kinase 7-1 (MPK7-1, KOG0660), calcium-dependent protein kinase 26-like (CDPK, KOG0032), and the probable protein phosphatase 2C 47 (PP2C, KOG0698), played a key role in this PPI map, which participated in protein phosphorylation and dephosphorylation under drought stress. The proteins in this map were also involved in signal transduction, water and ion transport, cell structure and division, and photosynthesis and energy metabolism.

Verification of identified phosphoproteins by Pro-Q diamond staining and tandem mass spectrometry.
To enhance the credibility of results, Pro-Q diamond staining and tandem mass spectrometry were used to validate the phosphoproteins identified in this study. The results from Pro-Q diamond and Coomassie brilliant blue (CBB) showed that 33 protein spots were phosphorylated at 6 and 24 h under drought stress (Fig. S4). MALDI-TOF/TOF-MS identification showed that these protein spots represented 24 phosphoproteins, further confirming the results from phosphoproteome analysis (Table S7). Four spots (1, 11, 14, and 24) belonged to three SCPL phosphoproteins: fructose-bisphosphatealdolase (chloroplastic, gi|357157399), ribulosebisphosphate carboxylase small chain PW9 (chloroplastic-like, gi|357146204), and ATPase beta subunit (chloroplast, gi|193075564). These phosphoproteins are related to photosynthesis and energy metabolism.

Discussion
Protein kinases/phosphatases participated in signal perception and transduction. The reversible phosphorylation modification event is regulated by protein kinase and phosphorylase, which could supervise multiple biological processes of organisms 34 . In this study, numerous protein kinases and phosphatases were SCPL, and played key roles in signal transduction (Figs 6 and S8). ABA is a crucial signal molecule in signal transduction responses to various abiotic stresses, including drought stress 35 . PKABA1 (gi|721630009) phosphorylated at Ser622 was identified in this study and found to be significantly up-regulated at 24 h (Table S2b). PKABA1, a member of the sucrose non-fermenting1-related protein kinase (SnRK) subfamily, acts as an intermediate in ABA suppression of GA-responsive gene transcription 36 . Some evidence suggested that over-expression of TaSnRK2S could enhance the tolerance of abiotic stressors, such as drought and salt stress, in A. thaliana 37 . In addition, protein phosphatase 2C (PP2C) was another key regulator and got involved in inhibiting SnRK's substrate phosphorylation which caused signaling cascade in the ABA signaling pathway 38 . PP2C-47 (gi|357147223) was found and down regulated at 6 h, which was in agreement with Sheard's study 38 : PP2C acted as a constitutive negative regulaor of SnRKs whose autophosphorylation is required for kinase activity towards downstream targets. Thus, it indicated that Bd21 responded positively to defense the drought stress with phosphorylation at 6 h (Supplemental Table S2a) In plants, the signals from the upstream elicitor receptors/sensors to the downstream mitogen-activated protein kinase (MAPK) substrates largely occur in three steps: the MAP kinase kinase kinase (MAPKKK), MAP kinase kinase (MAPKK), and MAPK cascade 39 . Some experiments have shown that various MAPK cascade members are involved in mediating signal transduction by transcription regulating P and N deprivation 40 . The activation of these MAP kinases is a sequential process through the phosphorylation of some conserved residues within an activation loop in their kinase domains 41 . In our study, three phosphorylated proteins were identified: MPK17 (gi|405778409), MPK7-1 (gi|405778403), and MPK3 (gi|405778399) ( Table S1). MPK17 had two phosphorylation sites (Ser465, Ser491) and the phosphorylation level at Ser465 had significant up-regulation at 24 h. A recent study showed that MAP kinases could serve as upstream activators of the stress marker genes under saline conditions 41 . A. thaliana MPK20 and proline dehydrogenase might interact with each other and be part of a functional hypoosmotic stress-signaling pathway 42 . Sequence alignment and conservation analyses demonstrated that the expression of MPK17 was high, similar to Arabidopsis thaliana ( Fig. 5; Table S5c). In addition, we identified some important CDPKs, such as CDPK26-like (gi|357125710; Table S2a) phosphorylated at Ser18, which supports the view that many of the autophosphorylation sites in CDPKs are localized at the N terminus.

Transcription factors involved in gene expression.
Under abiotic stress, many defense-related genes are regulated by transcription factors (TFs). In the present study, many TFs were phosphorylated (Table S2;  Table 2). For example, nuclear TFs Y subunit A-7-like isoforms X1 (NF-Y A, gi|721672572), NF-Y B-1 (gi|357122032), and NF-Y beta (gi|721661240) were phosphorylated at Ser202 (_GAS(ph)PANQTGNRE_), Ser4   (_PDS(ph)DNEDS GGGGGIGGGGNNK_), and Ser105 (_LDTLDS(ph)PK_), respectively. Both NF-Y B and NF-Y A proteins were up-regulated, whereas the NF-Y beta protein was down-regulated. Some studies involving maize, Arabidopsis, and wheat have proven that NF-Y B subunits may respond to drought stress 43,44 . Nelson et al. 43 reported tolerance to drought stress in transgenic maize, based on increasing ZmNFYB2 expression. The Bd21 phosphorylation levels were conserved ( Fig. 5; Table S5c). Thus, these TFs responded positively to abiotic stresses. Heat shock proteins (Hsps) and heat shock TFs (Hsfs) are involved in cellular responses to various sources of stress beyond heat stress 45 . The conformation of Hsps is affected by phosphorylation modification, such as the formation of active sites, flexibility, and inter-domain communication 46 . In other researches, HsfA1, as a main regulator, could enhance thermo tolerance and other stresses tolerance 47,48 . In our study, HsfA-2e-like (gi|357111341) was identified in Table S2a and HsfA-2e-like and Hsp81-1 (gi|357148345) have a complex relationship with other SCPL phosphoproteins (Figs 6 and S8). Hsfs recognized and combined with the Hsps gene promoter region of the heat shock element conservation motif. Thus, Hsfs may regulate Hsps transcription and gene expression to improve the stress tolerance of plants.

Phosphoproteins involved in water and ion transport.
In the face of various biotic and abiotic stresses, particularly drought stress, plants gradually adjust water and ion transport mechanisms. AQPs belong to a membrane protein family that rapidly transports water under osmotic pressure 49 . In our study, three AQPs were identified as belonging to the plasma membrane intrinsic protein (PIP) and nodulin 26-like intrinsic protein (NIP) subfamilies (Table S2a). All were up-regulated at 6 h under drought stress, which indicates that they had the same function in this process. NIP2-2 and PIP2-7 were phosphorylated at Ser271 and Ser282, respectively, in agreement with Guenther et al. 50 . Nodulin 26 phosphorylation at the C terminal tail (at Ser262) enhances water transport and is involved in the Ca 2+ signal pathway to accelerate regulation of the water channel switch. Thus, phosphorylation could influence the activity of these AQPs.
There are three types of H + -ATPase: F-type, V-type, and P-type. The P-type is a large family of molecular pumps that act as intermediates in ATP hydrolysis and ion transport or ion counter-transport 51,52 . Moreover, activation of the plasma membrane H + -ATPase is associated with phosphorylation of a penultimate Thr and binding to the phosphorylated C-terminus of 14-3-3 protein to control the stomatal aperture 53 . A recent study has showed that the stomatal aperture is a limiting factor in photosynthesis and plant growth, and the overexpression of H + -ATPase in guard cells is useful for promoting plant growth 54 . We also identified some proteins related to H + -ATPase. Plasma membrane ATPase (gi|357166497) and V-type proton ATPase subunit a1 (gi|357125740) were up-regulated at 6 h and phosphorylated at Thr881 and Ser688 (Table S2a), which enhanced enzyme activity and related to osmotic regulation. ROS detoxification and stress defense. Many biotic and abiotic stresses might induce plant cells to produce ROS, leading to the acceleration of lipid peroxidation and leaf damage 55,56 or beyond the antioxidant system 57 , particularly when subjected to drought stress. ROS scavenging or detoxification is a defense strategy in plants, involving a series of enzymes and proteins that reduce oxidative damage under drought stress 58,59 . In this work, one of the probable phospholipid hydroperoxide glutathione peroxidases (PHGPx, gi|357165189) was identified as an SCPL protein (Table S2b), which displayed up-regulated at 24 h. Phosphorylated PHGPx could become active to eliminate ROS and play a key role in the detoxification response to drought stress.
In addition, to sustain normal morphological structure and physiological function, plant cells must complete normal homeostasis, including electron transfer and oxidation reduction, and activate oxygen-scavenging to resist stress. In this study, we identified photosynthesis-and electron transport-related phosphorylated proteins that alleviate stress from drought. Chlorophyll a-b binding protein (chloroplast; LHCB, gi|357121229) was identified to be a specific SCPL protein at 6 h, and had significant up-regulation, as identified in salt stress (Table S2a;  Table 2). LHCB was regulated by the redox state of the plastoquinone mediating electron transfer between photosystems I and II, which involve reversible phosphorylation. Moreover, LHCB occurred redox dependent phosphorylation by specific thylakoid-bound kinases to affect the reconstruction of light harvesting complexes 60 . Furthermore, the changes in physiological and ultrastructure (Fig. 1) also implied that chloroplast structure was injured under drought stress and LHCB phosphorylation activity was weakened, which influenced photosynthesis and energy metabolism.
LEA proteins are made up of a group of dehydration-inducible protein families and their function is to help maintain the minimum cellular water requirement. Comparison of dehydrin DHN3-like (gi|357149091 , Table  S1) with a LEA protein (RAB-17) in maize 61 , showed the same sequence motif of -SSSSSSSEDD- (Table S1). Two LEAs were identified in SCPL proteins (LEA1, gi|721622359; LWAB19.3, gi|721639026) and up-regulated at Ser81 and Ser91. In maize cells, phosphorylation of LEA could facilitate binding to the targeted proteins and involve nuclear location 62 . Moreover, LEA is involved in ABA induction under drought stress.
In addition, molecular chaperones were involved in regulating plant growth processes and the chaperone complexes of Hsp90 could respond to ABA stimulation activitily 63 . Hsp81-1 (gi|357148345), a molecular chaperone member of the Hsp90 family, has the same motif with some identified sequence phosphopeptides, such as -TTEKEIS(ph)DDEDEEEK-in A. thaliana 64 , rice 65 , and wheat 63 (Table S5a). Hsp90 plays a key role in protein folding 66 and participates in signal-transduction networks and protein degradation and trafficking [67][68][69] . Moreover, in the present study, Hsp81-1 was up-regulated at 6 h ( Table S2a), suggesting that it could play an active role in response to drought stress.

A putative pathway of drought stress responses in Bd21 leaves through phosphorylation modification.
Based on the results from this study and previous studies, we propose a putative pathway of drought stress response in Bd21 seedling leaves through phosphorylation modification (Fig. 6). When plants were subjected to drought stress, membrane sensors sensed external environment changes, and transferred the chemical signals to the cell, which involved ABA, Ca 2+ , and MAPK signal pathways. Phosphate kinases/phosphatases were activated by the corresponding signal molecules. Then downstream TFs and translation factors were regulated by these phosphate kinases/phosphatases and they induced the production of drought stress response proteins such as LEA, Hsps, and LHCB. Massive accumulation of free radicals and ROS caused membrane lipid peroxidation, damaging the membrane systems, which caused the up-regulation of ROS scavenging-related phosphoproteins to defend against drought stress. At the same time, ABA and Ca 2+ concentrations triggered stomata closure to defend against water stress. The AQPs and transporters on membranes also regulated the balance of osmosis by phosphorylation/dephosphorylation.

Materials and Methods
Plant materials and drought treatment. Bd21 kindly provided by Dr. John Vogel et al. 24 at USDA-ARS, was used as plant material in this study. Bd21 seeds were surface sterilized by 15% sodium hypochlorite for 3 min, and rinsed 4 times in sterile distilled water. Seeds were submerged in water for 12 h at room temperature, and then transferred to wet filter paper to germinate at room temperature (22-25 °C) for 24 h. The uniformly germinated seeds were selected to grow in plastic pots containing Hoagland solution, in which was changed every two days. In a greenhouse, the experimental conditions included daily photo cycle of 16 h light/8 h dark (26 °C/18 °C) and 65-75% air humidity. Three biological replicates were performed under the same conditions. At the three leaf stage, the seedlings for different drought treatment times (0, 6, 12, 24 and 48 h) were cultivated in Hoagland solution containing 20% PEG6000 70 (ψ s = − 0.75) measured by VAPRO pressure osmometer (Wescor 5520, USA), whereas the seedlings for the control were only cultivated in Hoagland solution (ψ s = − 0.045). After treatment for different times, some sampled leaves for each biological replicates were used for physiological indicator measurement and ultrastructure observation, and the remaining leaves were kept frozen in − 80 °C for later protein extraction.
Seedling morphology, physiological parameter measurement and leaf ultrastructure observation. Phenotypic, physiological characterization and ultrastructure of Bd21 seedling leaves under the control and treatment conditions were assessed. RWC measurements were conducted according to Larbi & Mekliche 71 with minor modifications. Proline content and MDA content were measured, according to Bates et al. 72 and Zhao et al. 73 with minor modifications, respectively. Chlorophyll (chlorophyll a and b) were determined according to Zhang at al. 19 . Statistical analyses of the physiological measurements were conducted using an analysis of independent Student's T-tests by SPSS statistics software (17.0). For ultrastructural observation by TEM, conventional chemical fixation, dehydration and embedding of 2 mm × 2 mm leaves segments were performed, according to the procedure described by Thiel et al. 74 . For electron microscopic analysis with a H-7500 transmission electron microscope (HITACHI Company, Japan).
Scientific RepoRts | 6:35280 | DOI: 10.1038/srep35280 Protein extraction. Leaf total proteins were extracted according to the method of Lv et al. 21 with minor modifications. Approximately 500 mg of fresh leaves from each biological replicate was ground into a fine powder in liquid nitrogen. The ground powder was suspended in 4 mL SDS buffer (30% sucrose, 2% SDS, 100 mM Tris-HCl, pH 8.0, 50 mM EDTA-Na 2 , 20 mM DTT) and 4 mL phenol (Tris-buffered, pH 8.0) in a 10 mL tube, followed by the addition of 1 mM phenylmethanesulfonyl fluoride (PMSF) and 20 μ L protease inhibitor cocktail (Merck KGaA, Germany) per 1 g of fresh sample to inhibit the protease activities. Samples were mixed vigorously for 15 min at room temperature, centrifuged twice at 15700× g for 15 min each, and the supernatants were precipitated with 100 mM cold ethanolamine-methanol solution at − 20 °C over night. After centrifuging at 15700× g for 15 min, the pellets were rinsed with cold acetone (− 20 °C) and further centrifuged three times. After freeze-drying, the pellets were added to 300 μ L of solubilization buffer (7 M urea, 2 M Thiourea, 1% DTT (w/v) and 4% CHAPS) at room temperature for 2 h. Insoluble materials were removed by centrifugation at 15700× g for 15 min, and protein samples were measured using a 2-D Quant Kit (Amersham Bioscience, USA). The final protein solution was stored at − 80 °C for later use.
Phosphopeptide enrichment using TiO 2 microcolumns. Extracted protein mixtures were directly reduced with dithiothreitol (DTT), alkylated with iodoacetamide, and subsequently digested with endoproteinase Lys-C and trypsin as described previously 75 . The enrichment for phosphopeptides procedure was performed as reported by Wu et al. 76 with some modifications. In detail, the TiO 2 beads (GL Sciences, Tokyo, Japan) were incubated in 400 μ L loading buffer containing 65% Acetonitrile (ACN)/2% trifluoroacetic acid (TFA)/saturated by glutamic acid. A total of 2 mg tryptic peptide were resolved in 600 μ L loading buffer, and then incubated with appropriate amount of the TiO 2 beads. After washing by 600 μ L wash buffer (65% ACN/0.1% TFA), phosphopeptides were eluted twice with 300 μ L elution buffer (500 mM NH4OH/60% ACN). The eluates were dried down and reconstituted in 0.1% formic acid (FA)/H 2 O for MS analysis.
Phosphopeptide identification using LC-MS/MS. Enriched phosphopeptides were separated on a self-packed C18 reversed phase column (75 μ m I. D., 150 mm length) (Column Technology Inc.), which directly connected on the nano electrospray ion source on a LTQ-Orbitrap XL mass spectrometer (Thermo Fisher Scientific). Pump flow was split to achieve a flow rate at 1 μ L/min for sample loading and 300 nL/min for MS analysis. The mobile phases consisted of 0.1% FA (A) and 0.1% FA and 80% ACN (B). A five-step linear gradient of 5% to 30% B in 105 min, 35% to 90% B in 16 min, 90% B in 4 min, 90% to 2% B for 0.5 min and 2% B for 14.5 min was employed. The spray voltage was set to 2.0 kV and the temperature of heated capillary was 240 °C. For each sample, triplicate LC-MS/MS experiments were performed.
For data acquisition, each MS scan was acquired at a resolution of 60,000 (at 400 m/z) with lock mass option enabled and was followed by a data-dependent top 10 MS/MS scans using collision induced dissociation (CID). The threshold for precursor ion selection was 500. Mass window for precursor ion selection was 2.0 Da. The dynamic exclusion duration was 120 s, repeat count was 1 and repeat duration was 30 s. The analyzer for the MS scans was Orbitrap and for the MS/MS scans was LTQ (37% relative collision energy). Three biological replicates were performed from sample collection to the phosphopeptide identification using LC-MS/MS.
The raw files were processed with Maxquant (version 1.2.2.25) 77 and searched against own Brachypodium distachyon L. Database (NCBI_B. distachyon_25824 entries_20150320. Fasta). Up to two missing cleavage points were allowed. The precursor ion mass tolerances were 7 ppm, and fragment ion mass tolerance was 0.5 Da for MS/MS spectra. The false discovery rate (FDR) was set to < 1.0% for both peptide and protein identification, the minimum peptide length was set to 6.
Label free quantification and phosphorylation residue localization. Label free quantification was similar to the qualitative analysis on the experimental operation, but it needed more complex data requirements in the data extraction and processing methods. Each sample should be measured at least three times in parallel, in order to eliminate accidental errors. These datas of each phosphopeptide were normalized in calculated. MaxQuant software was used to calculated the integrating the ion intensities over its chromatographic elution profile coupled 26,77,78 . Phosphorylation residue localization was evaluated based on PTM scores that assign probabilities for each of the possible residues according to their residue-determining ions. In this study, Maxquant (version 1.2.2.25) was used to calculate PTM scores and PTM localization probabilities. Then, potential phosphorylation residues were grouped into three categories depending on their PTM localization probabilities 75 namely class I (localization probability, p ≥ 0.9), class II (0.9 > p ≥ 0.75), class III (0.75 > p ≥ 0.5) and class IV (p < 0.5). The FDR of 1% was used for phosphorylation residues identification. Spectra without residue-determining ions resulted in the identification of phosphopeptides with undetermined residues. Verification of identified phosphoproteins by Pro-Q diamond gel staining. Extracted leaf total proteins were separated by 2-DE and then 2D gels were stained with Pro-Q Diamond (Invitrogen, USA) according to the manufacturer's instructions. Briefly, the gels were fixed twice for 30 min/each time, and washed three times with ddH 2 O for 10 min/each time. The gels were incubated in Pro-Q Diamond stain in darkness for 2 h. The gels were destained with 20% acetonitrile in 50 mM sodium acetate (pH 4.0) for four times (30 min/each time), and then thrice washed for 5 min/each time. Finally, the gels were scanned on a TyphoonTM 9400 scanner (GE Healthcare, USA) with a 532 nm excitation laser and a 610 nm long pass filter with a resolution of 100 microns. After fluorescent image acquisition, the gels were stained with CBB to visualize total proteins. In addition, the same 2D gels were stained with CBB to visualize total proteins and used as control. The protein spots stained by Pro-Q Diamond were further identified by matrix-assisted laser desorption/ionization time-of-flight/ time-of-flight mass spectrometry (MALDI-TOF/TOF MS) according to Lv et al. 21 .
Scientific RepoRts | 6:35280 | DOI: 10.1038/srep35280 Bioinformatics. The significantly enriched phosphorylation motifs were extracted from the Motif-X algorithm 79 (http://motif-x.med.harvard.edu/). Search Tool for the Retrieval of Interacting Genes/Proteins (STRING version 9.1) database of physical and functional interactions was used to analyze the Protein-Protein Interaction (PPI) of all the phosphoproteins identified in the current study 80 (http://string-db.org). The sequences of all phosphoproteins were blasted in NCBI. Function annotations were based on Blast2Go 81,82 (http://www.blast2go.com/ b2ghome). Using EuKaryotic Orthologous Groups database (http://eggnog.embl.de/version_3.0/), KOG search and numbers were obtained. Then a dataset containing all the KOG numbers was used to do PPI analysis by using STRING and the interested network was displayed by Cytoscape (version 3.1.0) software (http://www.cytoscape. org/). Secondary structure of phosphoproteins was predicted via SABLE server version 2 (http://sable.cchmc. org/). Phyre2 online server (/redict the 3D structure of some interested phosphoproteins. Then the phosphorylated residues were displayed by SPDBV (version 4.1) software 83 (http://spdbv.vital-it.ch/).