Proteomic profiles of cotton fiber developmental 2 transition from cell elongation to secondary wall 3 deposition 4

Cotton fiber development transition from elongation to secondary cell wall biosynthesis 15 is a critical growth shifting phase that affects cotton fiber final length, strength and other properties. 16 Morphological dynamic analysis indicates that an asynchronous fiber developmental pattern 17 between two cotton species. The critical time point for Gh and Gb fiber elongation termination is, 18 respectively, 23 and 27 days post-anthesis (dpa). The temporal changes of protein expression at three 19 representative development periods (15-19, 19-23, 23-27 dpa) were examined in both species with 20 iTRAQ technics. Strikingly, a large proportion of differentially expressed proteins (DEPs) was 21 identified at 19-23 dpa in Gh or at 23-27 dpa in Gb, corresponding to their fiber developmental 22 transition timing from elongation to secondary cell wall biosynthesis. To better understand fibers 23 transitional development, we comparatively analyzed those DEPs in 19-23 dpa of Gh vs. in 23-27 24 dpa of Gb, and noted that these cotton species indeed share fundamentally similar fiber 25 development features under the biological processes. It also showed that there have limited overlaps 26 in both specific upregulated and downregulated proteins between the two species, suggesting 27 specie-specific protein regulations in development. Proteomic profiling revealed dynamic changes 28 of several key proteins and biological processes that potentially correlate with fiber development 29 transition. During the transition, upregulated proteins mainly involved in carbohydrate/energy 30 metabolism, oxidation-reduction, cytoskeleton, protein turnover, Ca2+ signaling etc, whereas 31 important downregulated proteins mostly concentrated in phenylpropanoid and flavonoid 32 secondary metabolism pathways. Several changed proteins in this key stage were also validated by 33 qRT-PCR. Overall, the present study provides accurate pictures of the regulatory networks of 34 functional proteins during the fiber developmental transition. 35


Introduction
Cotton (Gossypium spp.) is one of the most important crops, contributing prevalent natural textile fiber worldwide.The most commonly cultivated cotton species today are upland cotton (G.hirsutum.,Gh) and sea island cotton (G.barbadense, Gb), both of which are tetraploid plants originated from interspecific hybridization event about 1-2 million years ago [1].Selective breeding of the cottons has led to the high yield and diverse environmental adaptability of Gh, therefore the Gh cultivation now accounts for the majority of cotton fiber production, whereas the Gb grown in selected environments is prized for superior fiber length, strength, and fineness [1].Nowadays, within the premium textiles market, there is a demand for higher-quality, high-yield cotton fibers, however the molecular mechanics governing fiber quality are still not well understood.Unraveling the molecular basis for different fiber agronomic traits between Gh and Gb will contribute significantly to the characterization and manipulation of the specific genes that control fiber quality and yield, thereby allowing for the improvement of cotton fibers.
Cotton fibers are highly elongated and thickened seed epidermis single cells that undergo four major sequential and overlapping developmental stages: fiber initiation, elongation (primary cell wall synthesis), cell wall thickening (secondary cell wall deposition) and maturation [2].Given that the cotton fiber quality is determined by the final length and strength, many endeavors have devoted to investigating the regulatory mechanisms underlying the fiber cell elongation and secondary cell wall deposition during cotton development.In recent years, knowledge of cotton genome structures, interrelationships between cotton varieties in relation to development and evolution have obtained impressive increases by using multi-omics analysis approaches [3].
Proteins are the direct performers for most biological activities and functions.Proteomic analysis that provides overall information about protein regulation and active pathways has been widely used to reveal molecular mechanisms of particular biological processes [4].Amongst the various proteomic technologies, the "isobaric tags for relative and absolute quantitation" (iTRAQ) method can identify numerous proteins and provide more reliable quantitative information than conventional analysis by two-dimensional gel electrophoresis [4].Prior studies using iTRAQ analysis have gained insight into the differences in gene/protein expression of cotton in response to environment stress [5], and domestication [6], as well as fiber development [7][8][9].
The fiber developmental transition stage is often thought of as a period when primary and secondary wall deposition overlap [10].During the transition, fiber development is accompanied with significant changes in physiological processes and cell wall protein contents, requiring organization and rearrangements of polysaccharides.Here, we have identified a significant variation at the timing of development transition among different cotton species by fiber morphological and proteomic dynamic analysis.Based on comparative proteomic analysis, we revealed two cotton species, Gh and Gb, share a highly similar development regulatory patterns during their respective fiber cell transition.Furthermore, several key interspecific differentially regulated proteins that are potentially involved in the cotton fiber developmental transition stage were identified and analyzed.
This study provides new clues concerning the fiber development transition at proteomic level, thereby highlighting candidate genes/proteins and related pathways for cotton fiber improvement.

Measurement of fiber length and thickness
Two cotton cultivars, Gh cv.Xinluzao 36 and Gb cv.Xinhai 2, were healthily grown in fields at the Xinjiang Horticulture Experimental Station.Flowers were tagged at anthesis, and developing bolls were harvested every 4 days in an interval from 11 to 35 dpa.The bolls were dissected immediately, and the length of fibers were determined using a previously reported method [11].To determine the fiber thickness, fibers were fixed in 3% glutaraldehyde and dehydrated in an ethanol series (from 30% to 100%) before then being infiltrated with Spurrs resin (Electron Microscopy Sciences).The thickness of the cell wall was examined by measuring the cross-section fibers under a transmission electron microscope (TEM) (Bio-TEM H-600, Hitachi, Japan).All measurements were conducted with at least 100 fibers and 10 different ovules at corresponding boll age.

Protein extraction of cotton fiber
Cotton fiber protein extractions were performed using a modified phenol extraction method as reported [12].For protein extraction, a total of 24 independent Gh and Gb samples were collected at time points 15, 19, 25, and 27 dpa.Approximately 800 mg fibers were finely ground with liquid nitrogen, and then was homogenized in 5 ml buffer (50 mM Tris-HCl, pH 8.0, 30% sucrose, 2% SDS, 1% DTT, and 1mM PMSF).After adding saturated phenol, the mixture was vortexed and centrifuged at 10,000 ×g for 10 min at 4 ℃.The upper phenol phase was collected and mixed with 5 ml NH4AC (0.1 M) at -20 ℃ for 30 min, and then centrifuged at 20,000 ×g for 20 min at 4 ℃.The pellets were collected and washed with cold NH4AC (0.1 M), and 80% acetone to obtain proteins.The proteins were vacuum dried and stored at -80℃.

Proteome analysis of cotton fibers by iTRAQ
Fiber proteomic data were collected at the National Center for Protein Sciences Beijing, China.
Before labeling, the total protein of each sample (100 μg) was digested and reconstituted using 8-plex iTRAQ reagent (AB Sciex Inc., CA, USA).The fiber proteins of Gh at the developmental time point of 15, 19, 23 and 27 dpa were labeled with iTRAQ tags 114, 116, 118 and 121, and those of Gb were labeled with iTRAQ tags 113, 115, 117, and 119, respectively.Chromatography consisted of Thermo Surveyor HPLC system was operated at 500 nL per minute via a split solvent line.Each sample was loaded on a BioBasic C18 reversed phase column (Thermo 72,105-100,266) and flushed for 20 min with 5% acetonitrile (ACN), 0.1% formic acid to remove salts.Peptide separation was achieved using a Thermo Surveyor MS pump with a gradient HPLC method washing from 5% ACN to 50% ACN in 620 min, followed by a 20 min wash of 95% ACN and equilibration with 5% ACN for 15 min.The Surveyor was coupled with a Thermo LCQ DECA XP Plus mass spectrometer with a stock nanospray ion source.Data were acquired with a 2.5 kV ion spray voltage, 30 PSI curtain gas, 5 PSI nebulizer gas, and 150 ℃ interface heater temperature.Each cycle time was fixed to 2.5 s.Dynamic mass exclusion windows were 2 min long with a repeat count.All samples were run with three replicates.

Protein identification and bioinformatical analysis
Protein identification and quantification were performed using the Mascot 2.3.02software (Matrix Science, Boston, USA).MS/MS spectra were analyzed with Protein Pilot software (Protein Pilot 4.0; AB SCIEX) against the corresponding genome databases using the Paragon algorithm.An automatic decoy database search strategy was employed to estimate the false discovery rate (FDR) using the PSPEP software integrated in the Protein Pilot.Proteins were identified using the following parameters: sample type = iTRAQ 8-plex (peptide-labeled), Cys; alkylation = iodoacetamide; digestion = trypsin; instrument = Triple TOF5600 (AB SCIEX).To annotate coding sequences with the highest score, we searched the non-redundant protein sequence database at NCBI.Differential expression from the protein data was judged with the following criteria: number of unique peptides≥2; threshold of fold change for upregulation/downregulation = 1.5/0.67;and maximum allowed fold change = 30.The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis were performed with identified differentially expressed proteins [13,14].Once p < 0.05, the GO term or pathway was regarded as a significant enrichment.

Quantitative real-time PCR (qRT-PCR)
Total RNA was isolated from fiber samples using a Total RNA Isolation Kit (Biorbyt, San Francisco, United States).First strand cDNAs were synthesized from 1 μg of total RNA using the ReverTra Ace qPCR RT kit (TOYOBO, Osaka, Japan).The reverse transcription product was diluted 30-fold with RNase-free water and stored at -80 °C.The specific primers for selected genes were designed and synthesized by Sangon Biotech Co., Ltd.(Shanghai, China).Cotton Ubiquitin7 gene was used as a reference gene to normalize the cDNA amplification in each reaction.Triplicate replicates of qRT-PCR were performed with SYBR Premix Ex Taq (TaKaRa, Dalin, China) on ABI 7500-Fast Real Time PCR system (Applied Biosystems, CA, USA).Relative gene expression levels were calculated using the 2 -ΔCT method.The amplification primers are listed in Supplemental Table 2 Cotton fiber cells undergo substantial elongation and expansion throughout development.To investigate the interspecific divergence of fiber development between Gh and Gb, we determined the fiber lengths across from 5 to 37 dpa.The fiber length approximately linearly increased over the 5-19 dpa in Gh or 5-23 dpa in Gb, respectively.The elongation of Gh fiber almost ceased after 23 dpa, while the elongation in Gb fibers ceased around 27 dpa (Fig. 1A).The elongation rate (length increase per day) within each species declined sharply after 19 dpa, whereas a higher fiber elongation rate occurred in Gb at later periods (Figure 1B).This identified prolonged fiber elongation in Gb cotton, consistent with prior reports [10].Plant cell wall thickening usually is the effect of secondary cell wall biosynthesis.Fiber cell wall thickness was also examined using TEM.It showed that the fiber cell wall thickness has no apparent difference between two cotton species during the fibers rapid elongation period.Fiber cell wall thickening initiation in Gh cotton was at 19 dpa, beyond that time point the cell wall thickness dramatically increased, but the onset of cell wall thickening in Gb fibers had a delay, starting at around 23 dpa (Figure 1C).Compared with Gh, a faster thickening rate was also found in Gb during subsequent fiber growth ( Fig. 1D).

Morphological dynamic analysis of cotton fiber development
These results suggested that two respective periods, 19-23 dpa and 23-27 dpa, were crucial for Gh and Gb development respectively, involving the transition phase from elongation to secondary cell wall deposition.Prior studies of the two near-isogenic cotton lines showed that the transition to secondary wall deposition correlates with their fiber bundle strength differences, and the duration of this transition stage may determine cotton fiber length, as well as other properties [15].Importantly, the prominent phenotypic differences between these two cotton species does suggest an asynchronous fiber developmental pattern in transition from elongation to secondary cell wall deposition.

Proteomic analysis of cotton fiber development
To shed light on the mechanisms controlling fiber development, we examined proteome changes during fiber rapid elongation and structurally thickening stages.The eight-plex iTRAQ experiments allowed for a detailed comparison of cotton fiber protein expression differences at adjacent periods (15-19, 19-23, and 23-27 dpa).A total of 1197 proteins were successfully identified at a 95% confidence level and a 1.0% FDR.Finally, we quantified and annotated 797 distinct proteins with two or more unique peptides by using cotton databases [16,17].Proteins that showed a difference in abundance corresponding to at least a 1.5-fold change and a P value of < 0.05 were considered to be differentially expressed proteins (DEPs).Based on these criteria, 112 and 94 DEPs were identified in Gh and Gb, respectively.A total of 102 unique DEPs (12.7 % of 797 proteins) among two cotton species were identified at those adjacent developing time points To classify these DEPs, we performed the GO term enrichment analysis, which was divided into molecular function, biological processes, and cell composition (Fig. 2 and Supporting Information in previous studies on cotton fiber development [9].In biological process category, organic substance metabolic process (59.7%), primary metabolic process (56.8%), cellular metabolic process (55.5%), single-organism metabolic process (41.6%) presented the most significant four enrichment among others.Given that the enriched proteins mainly involved in various metabolic processes, the regulation of the basal metabolic reactions is believed to play a critical role in sustaining a rapid developmental growth of cotton fiber cell.In the cellular component category, cell part (52.3%) was the most abundant subcategory, followed by intracellular (50.7%), membrane-bounded organelle  To comprehend pathways in the fiber development, we also preformed KEGG enrichment analysis on those DEPs.It showed that several key metabolism and protein-related biochemical pathways were significantly enriched (FDR < 0.01) (Fig. 3).The enrichments are mainly associated with metabolic pathways (30.8%), biosynthesis of secondary metabolites (22.1%), carbon metabolism (13.6%), ribosome (12.0%), glycolysis/gluconeogenesis (8.4%), citrate cycle (TCA cycle) (6.2%), glutathione metabolism (6.2%), biosynthesis of amino acids (6.2%), pyruvate metabolism (5.8%), protein processing in endoplasmic reticulum (5.2%), proteasome (4.2%) and peroxisome (3.6%) (Fig

3.
).The distribution of the differentially regulated proteins indicates a strong role of energy/carbohydrate metabolism pathway throughout fiber development, as well as remarkable secondary metabolic pathway regulations involved.This result is plausible, because rapid cell elongation and fiber development require a large amount of energy and carbon intermediates for cell wall synthesis [18].During this distinct transition stage, cotton fiber development experiences some important physiological and biochemical events occurring, for example, changes of metabolic sugar contents [20], the degradation of the cotton fiber middle lamella [21], the deposition of winding cell wall layer [22].Therefore, it is presumed that there must be a marked change in protein expression pattern associated with the transition stage of fiber development [19].Identification of the significant proteomic variations within cotton fibers switching to secondary cell wall deposition indicates that the critical transition stage can be distinguished by protein expression dynamics.Based on these analyses, Gh and Gb cotton undergo the asynchronous fiber developmental transition process at protein level.Cotton fiber development transition is a significant shifting phase and is worthy of further investigation to discover critical developmental factors responsible for fiber length, stiffness and strength [23].

Comparative analysis of differentially expressed proteins during fiber developmental transition between two cotton species
The fiber development transition experiences extreme morphological changes from elongation to cell wall thickening and generates significant protein expression variation among Gh and Gb.To gain insight into the cotton fiber development transition, the DEPs from these distinct periods for both species (19-23 dpa for Gh vs. 23-27 dpa for Gb) were profiled and compared on basis of biological process.In the stage, 31 upregulated and 26 downregulated proteins were identified in Gh cotton, while 27 upregulated and 21 downregulated members were recorded for Gh cotton fibers.Notably, we found that Gb and Gh indeed shared a highly similar development when these two stages, 19-23 dpa for Gh and 23-27 dpa for Gb, were compared under the biological processes.Function comparative analysis revealed changes in several key biological processes, identified by DEP number in each category, with significant similarity between the two species.In both cottons, the upregulated proteins were mainly involved in carbohydrate metabolism, oxidation-reduction, cytoskeleton organization, response to calcium ion, proteolysis, glycolysis, lipid metabolism, signal transfection, protein folding and transport process, whereas the downregulated proteins were mainly involved in ribosome biogenesis, secondary metabolism, signal transduction, nucleic acids processing, protein folding and amino acid metabolic processes, among others (Fig. 5).Furthermore, a pairwise comparison with DEPs that identified in the two cotton species respective transitions was undertaken to explore protein expression patterns.All identified DEPs from this key development stage were summarized and listed in Tables 1. Eighteen common upregulated and 15 common down-regulated proteins were found in two cotton species (Fig. 6).These common proteins account for 58.1%, 66.7% of upregulated proteins and 57.7%, 71.4% of downregulated proteins in Gh and Gb, respectively.This suggests that two cotton species have only partial overlap in their fiber-developmental proteome during the transition, with some protein regulation being unique to both species during fiber development.This expressed protein variation could be responsible for interspecific phenotypic differences that include, for example, fiber traits.

Key gene/proteins and pathways involved in the cotton fiber developmental transition
During the developmental transition between primary and secondary wall deposition, the cotton fibers undergo spatial and temporal cell wall remodeling.In our identified protein categories, the largest number of upregulated proteins were involved in carbohydrate metabolisms.
Carbohydrate metabolism pathway provides the essential carbon skeletons for the synthesis of cell wall polysaccharides and fatty acids, as well as energy storages [18].Several remarkably upregulated common proteins in carbohydrate metabolism pathways have been identified in both Gh and Gb cotton species, including sucrose synthase (SUS), sucrose synthase-like (SUSL), cellulose synthase 8 (CESA8), endoglucanase (EG), and pectinesterase (PME), which have all been reported to be involved in cell wall biogenesis and important for cotton fiber production and quality [24,25].It is noting worth that of those DEPs, the PME is one of the most substantially regulated proteins, having a 11.29 and 17.45-folds increase in Gh and Gb, respectively (Table 1).Indeed, plant PME catalyzes the deesterification of pectin, one major component of the primary cell wall and middle lamella, has been shown to play crucial roles in regulating cell wall expansion, elongation and adhesion [25].Compared with Gh, a higher-fold PME upregulation occurred in Gb fibers at the transition phase, consistent with prior reports that a higher PME activity was found during the later stages of cotton fiber development [25].In Gb cotton fibers, other several key carbohydrate metabolism pathway enzymes, SUS, SUSL, EG and CESA8 also showed much more folds expression upregulation, suggesting Gb fiber transition has a superior carbohydrate metabolism than the transition in Gh.Meanwhile, several proteins distinct to both species were detected in proteomic analysis, such as upregulated UDP-D-glucose dehydrogenase in Gb, and reinforced expressed acid beta-fructofuranosidase-like protein in Gh.
These results suggest that stage-specific, upregulated critical enzymes in carbohydrate metabolism pathways would facilitate fiber cell wall developmental transition and might contribute to the variation in fiber traits between species.Furthermore, secondary wall cellulose micro-fibril formation in cotton fiber cells is an energetically costly process [2].Multiple proteins involved in glycolysis, TCA cycle, including enolase, ADP/ATP carrier proteins, pyruvate dehydrogenase E1 component subunit beta-3 and glyceraldehyde-3-phosphate dehydrogenase, presented in the up-regulation category, suggesting that 'energy' production is still active cellular process at fiber transition stage and might be a basis of the fiber physiological changes.The N.D represents that not obvious protein abundance difference was detected in the assays.
Differentially expressed proteins are involved in oxidation-reduction processes during cotton fiber development [10].During the fiber developmental transition, ascorbate peroxidase (APX) family members, APX and APX6, and NADP-dependent D-sorbitol-6-phosphate dehydrogenase-like protein (S6PDH) were found to be markedly upregulated in both cotton species.Alteration of fiber cell reactive oxygen H2O2 levels in in-vitro ovule cultures has been reported to affect the differentiation of the cotton fiber cell wall [26].The increase in reactive oxygen (ROX) scavenging enzymes APX and APX6 could maintain a low H2O2 level and regulate intracellular reactive oxygen species homeostasis, thus indicating the importance in regulating H2O2 related signal pathways for cotton fiber development Alternatively, the fiber cell elongation protein Ghfe1, as oxidationreduction related protein [27], appeared downregulated in both cottons, is consistence with the fiber elongation gradually ceased during the transition.We also noted Gh cotton has upregulated the peroxygenase 3 and the monodehydroascorbate reductase-like protein, as well as Gb specie has Protein turnover is the net result of continuous synthesis and breakdown of body proteins and ensures maintenance of optimally functioning proteins in organisms [28].Protein synthesis, folding and degradation pathways therefore are associated with protein turnover and amino acid biosynthesis.When the major mass of the fiber becomes crystalline cellulose, the total protein content of the developing cotton fiber eventually decreases during the fiber development process [10].
Among these identified proteins, the expression of 'protein degradation' class, including aspartic protease-1 like, proteasome subunit beta type and carboxypetidase Y-like protein, are increased in the development transition.Meanwhile, the 'protein synthesis' class that consisted of largely of ribosomal proteins and the heat shock protein family were down-regulated in this stage, such as these 60S ribosomal protein L8-3-like; 40S ribosomal protein S14-3-like and 60S ribosomal protein L18a, and 20 kDa chaperonin-like underwent obvious downregulation in both species.This turnover regulation would allow plants to balance protein synthesis and degradation during fiber developmental transition.
Fiber cell morphology is largely determined by the highly dynamic cytoskeleton architecture.In this distinct stage, several cytoskeletal-related proteins with a dramatic increase, including the βtubulin-3 (TUB-3), annexin (ANN) were found in both cotton species.Upregulation of cytoskeletalrelated proteins are essential for the fibers normal morphogenesis changes.Indeed, these ANN and TUB-3 have been reported to participate in cytoskeleton dynamic assembling and maintenance [29,30].Interestingly, other two known actin dynamic regulators, profilin and acid cyclase-associated protein (CAP) are differentially upregulated among two cotton species.Dynamically changed profilin that binds actin cytoskeleton has been shown to be important for restructuring cell shape [31].However, we detected the increased profilin only in Gh but not in Gb fibers.Additionally, a greater than 2-fold upregulation of another actin-binding protein CAP was found in Gb cotton.
Identification of these cytoskeletal related proteins in development transition further supported reorganization of the actin cytoskeleton is an important scheme in controlling direction of cellulose fibril deposition in the developmental process switching from elongation to secondary wall deposition.
In particular, among these DEPs, it appears that Ca 2+ -signaling pathways are involved in fiber development transition.Calcium mediated signaling plays an important role in cell division and differentiation including root hair elongation [32].Preferential expression of calcium binding proteins during fiber initiation and elongation stages have been reported in cotton [10].Several highly upregulated Ca 2+-signaling pathway proteins, including calmodulin-7 and calreticulin-like, were detected in the developmental transition for both cottons.A few recent reports have shown the overexpression of Ca 2+ -dependent protein kinase1 (CPK1) stimulates the onset of secondary wall deposition [33].The probable calcium-binding protein CML13 was preferentially upregulated in Gh fibers, which might regulate Ca 2+ homeostasis in the developing fiber cells.
Several other signaling transduction molecules were identified in cotton fiber development transition.Among them, ERBB-3 binding protein 1 and 14-3-3-like protein 2 were commonly downregulated in both species.It has reported that 14-3-3 proteins participated in regulation of fiber initiation and elongation by modulating brassinosteroid signaling and overexpression promotes fiber elongation in cotton [34].Therefore, the downregulated 14-3-3 protein is plausible actor in the developmental transition from fiber elongation to secondary cell wall synthesis.
Additionally, in the fiber developmental transition stage, 14-3-3 protein 6-like, fasciclin-like arabinogalactan protein 1 and STS14 protein-like were only downregulated in Gh, whereas another key signal transduction related molecule, Ras-related protein RABA1f-like, was downregulated in Gb.Moreover, a sulfurtransferase is reported to be abundant in the plant cell wall and plasma membrane, promoting cell differentiation [35].In Gh cotton fibers, the upregulation of sulfurtransferase as well as guanosine nucleotide diphosphate dissociation inhibitor (GDI) was recorded.While another upregulated signal transduction molecule frasciclin-like arabinogalactan protein 1 (FLAs1) was found only in Gb.These results indicate that several species-specific signal transductions may take place during this process and altered expression of various signaling molecules underpin a complicated signal transduction network regulating many downstream biochemical events during the transition phase.
While cellulose synthesis is gradually prevails during the late stages of fiber development, many metabolic pathways that are active during fiber elongation are repressed.In this developmental transition, multiple secondary metabolism related proteins, including shikimate Ohydroxycinnamoyltransferase-like (HCTL), naringenin (NAR), anthocyanidin reductase-like (ANRL) and chalcone-flavanone isomerase family protein (CHI), were found to be significantly downregulated in both species.Phenylpropanoids are a group of plant secondary metabolites derived from phenylalanine which have a wide variety of functions both as structural and signaling molecules [36].Plant HCTL enzyme participates in phenylpropanoid biosynthesis [37].The marked downregulation in HCTL expression, consistent with prior reports that cotton fiber thickening is negatively related to phenylpropanoid content, suggest that the HTCL associated with phenylalanine metabolic pathways might play a key role in the fiber developmental transition process.Also representing the temporal shift in protein expression pattern are down-regulated flavonoid biosynthetic pathways.As Yoo et al. reported, carbon resources might be reallocated in developing cotton fibers [3].In this proteomic analysis, three enzymes, NAR, ANRL and CHI, all involved in the flavonoid biosynthesis pathway [38], were downregulated throughout this development stage in both species, suggesting that carbon resources have been transferred away from flavonoid metabolism and might be coordinated to direct carbon flux into cellulose during this stage, thus highlighting that phenylpropanoid and flavonoid metabolism represents a novel pathway with potential for cotton fiber improvement.

Validation of differently expressed proteins by qRT-PCR
To examine whether the differences in abundance were consistent with the differences at the mRNA levels, qRT-PCR was used to analyze the transcripts encoding 12 proteins with differing expression, belonging to various functional categories, at least two time points during fiber development stages (15-19, 19-23, and 23-27 dpa).
The qRT-PCR analysis indicated that expression patterns of 9 mRNA abundance (75% of 12 proteins) were highly consistent with the protein data, while the remained 3 proteins (25%) were partially consistent (Fig. 7).Nevertheless, the transcript and protein levels of ANN and CHI conflicted at 19-23 dpa.As has been noted previously, mRNA levels do not always correlate well with the level of corresponding protein mainly because of post-transcriptional regulation mechanisms such as nuclear export and mRNA localization, transcript stability, translational regulation, and protein degradation [39].These results confirmed the reliability of the proteomic analysis results and proteomic analysis is essential for identifying the final products responsible for different cellular functions.Proteomic analyses are crucial for providing accurate pictures of the regulatory networks of functional genes/proteins.

Figure 1 .
Figure 1.The dynamic change of fiber length and cell wall thickness during development.(A).Gh and Gb cotton fiber length at different dpa.Standard deviation (SD) of error bars were calculated with fifteen biological replicates.(B).Gh and Gb cotton fiber elongation rate across development time.(C).Gh and Gb cotton fiber cell wall thickness at different dpa.Approximately 30 fibers were measured for each sample.(D).Gh and Gb cotton fiber cell wall thickness rate across time.

Figure 2 .
Figure 2. Gene Ontology analyses of differentially expressed proteins during fiber cell development.The ten most predominant enriched terms are shown in the GO categories.

PreprintsFigure 3 .
Figure 3. KEGG pathway enrichment analysis of differentially expressed proteins during the fiber cell development.The twelve most significantly enriched pathways were shown.
Upon further inspection of the identified DEPs across 15-19, 19-23 and 23-27 dpa for both species, we noted that the DEPs quantities are unequally distributed between the adjacent time points in each cotton specie (Fig. 4).Twenty (15-19 dpa), 57 (19-23 dpa), and 35 (23-27 dpa) DEPs were identified in Gh, and 23 (15-19 dpa), 23 (19-23 dpa) and 48 (23-27 dpa) DEPs were found in Gb, respectively.There is a distinct developmental period in both cotton species where the identified DEPs were almost double in comparison with those of all other periods, representing a burst stage in protein differential regulation (Fig 4.).Interestingly, the identified largest amount of DEPs was between 19 and 23 dpa in Gh, having a total of 57 DEPs, whereas the maximum number of DEPs was instead of presenting at the same time intervals, and that happened between 23 and 27 dpa in Gb, with 48 DEPs (Fig 4.) Strikingly, the delayed burst in DEPs in Gb is consistent with the morphological dynamic analysis showing that Gb fiber development exhibited a delayed onset in secondary wall deposition.Moreover, it has been reported that there have a large number of cotton gene regulations at transcription level in fiber developmental transition, usually in 19-27 dpa [19].

PreprintsFigure 4 .
Figure 4. Temporal changes of protein expression among two cotton species developing fibers.Numbers that beside of the arrow line designate the number of up-regulated/down-regulated proteins (at least 1.5-fold and FDR < 0.05) relative to their adjacent developmental time points.The highest number of differentially expressed proteins occurred in the transition to secondary cell wall synthesis.

Figure 5 .
Figure 5. Comparative analysis of differently expressed proteins under biological processes among respective protein regulation burst periods of two cotton species.The identified differently expressed proteins were classified into up-regulated and downregulated protein sets, and then carried out the GO biological processes annotation.

Figure 6 .
Figure 6.Venn diagram analysis of the differentially expressed proteins between Gh and Gb cotton during fiber development transition from elongation to second cell wall deposition.Venn diagrams showing the number of differentially expressed proteins and the overlap of identified common regulated proteins among Gh and Gb cotton during the fiber cell wall developmental transition.

Table S1 )
. Biological terms having significantly enriched (false discovery rate, FDR-corrected < 0.01) were shown.In molecular function class, ion binding (36.0% of the total number of DEPs) and small