- Research
- Open access
- Published:
Transcriptome dynamics along with expression analysis of key genes involved in fiber development between Gossypium barbadense and Gossypium darwinii
BMC Plant Biology volume 25, Article number: 691 (2025)
Abstract
The transcriptome profiling for underpinning the role of key genes controlling formation of fiber in cultivated Gossypium barbadense compared to wild allotetraploid cotton Gossypium darwinii which remained less investigated. Owing to excellent fiber quality of both Gossypium barbadense and Gossypium darwinii and information obtained via Simple Sequence Repeat (SSR) markers, two lines: Xh-18 and darwinii 5–7 were selected for transcriptome sequencing during developmental stages i.e., fiber initiation, elongation, and secondary cell wall (SCW) biosynthesis followed by 0 day after anthesis (DPA), 5DPA, 10DPA, 15DPA and 25DPA, respectively. Twelve libraries of RNA-seq were generated and sequenced individually generating approximately 818 million clean reads of Gossypium darwinii. However, for Gossypium barbadense more than 844 million clean reads were recorded. The Pearson Correlation Coefficient (PCC) analysis results indicated that gene expressions for both Gossypium barbadense and Gossypium darwinii indicated more than 90% of the commonalities at the same stage of fiber growth. However, genes found among Gossypium darwinii at 5, 6 and 7 DPA and XH-18 at 10 and 25 DPA were found disimilar. The expression quantity of RNA sequencing data, 31 genes were found common throughout all stages of DPAs in Gossypium darwinii 5–7 whereas, 377 genes were common in Gossypium barbadense XH-18 at 0, 5, 10, 15, 20, and 25 DPA stages of fiber development. Three genes XLOC_080616 (LTPG2_ARATH/ NLTL2_ARATH; GPI-anchored 2 non-specific lipid transfer protein, similar to At2 g13820), XLOC_065471 (LPAT2_ARATH/LPAT2_BRAOL 1-acyl-sn-glycerol-3-phosphate acyltransferase2), and XLOC_077416 (uncharacterized protein) shown up-regulated expression during 15 and 20 DPA for both Gossypium barbadense and Gossypium darwinii. It will also explore the possible role Gossypium barbadense DNA segment associated with development of fiber quality. Furthermore, this research will decipher the underlying process of fiber development and the possible role of genes for fiber formation in both Sea Island and wild cotton species.
Introduction
In Cotton, the fibers are highly specialized single-celled, unbranched and simple trichomes which arise from approximately 25% of the epidermal cells from the external covering of an emerging seed. It has been observed in Arabidopsis thaliana that various transcription factor and differentially expressed genes (DGEs) regulates the development of trichomes [6, 24, 46]. Fiber development in cotton has four partially coinciding stages namely, initiation, elongation'primary cell wall formation and maturation [13, 14]. The elongation of fiber cells take place between 10 and 16 DPA and lasts until 20 DPA. The primary and secondary cell wall synthesis (SCW) stage coincides within the duration of 16 to 25 DPA. The cell elongation speed declines and even ceases during the SCW formation stage. The production of cellulose plays an important role in fiber cells because it accounts for more than 95 percent of the dry weight of mature cotton fiber. Cotton fiber may be used as an outstanding model system for studying the cell development in plants [48]. The fiber cells length ranges from 30–40 mm whereas, approximately 15 mm range in thickness among cultivated cotton [26, 28, 48]. Cotton fiber initiation is an intricate biological procedure including fiber cell protrusion from ovule surface and early expansion, while the knowledge about regulatory mechanisms is vaguely understood. The lint fiber develop before fuzz fibers but initiation of fuzz fibers development timing varies among genotypes [37, 38]. During the growth of a fiber, the cells undergo a 1000 to 3000-fold elongation after initiation. The elongation period involves the increase in growth stage of the fiber, which might prolongs form 25 to 30 DPA [40].
During elongation, fibers form a thin, extensible primary cell wall consisting of several distinct carbohydrate polymers [33]. Normally, fiber elongation stage lasts from 5 to 25 DPA at maximum growth rate of more than 2 mm per day at 10- 12 DPA. Before transitioning to secondary cell wall production on 16 DPA, each cotton fiber stretches from 10 to 15 µm to 2.5 to 3.0 cm [35]. Cotton fibers enter the stage of SCW synthesis at 15 DPA, to produce a substantial amount of cellulose and thicker cell walls. Furthermore, cellulose production is so intense during the secondary cell wall development stage that cellulose accounts for approximately 95 percent of the dry weight in mature cotton fiber. Consequently, fiber traits i.e., micronaire, strength and length are determined by fiber cell elongation, process in which cell wall expands [16, 17, 24]. The tremendous development of high-throughput sequencing technology is blazing new paths for deciphering the mechanisms behind complex traits.
The whole genome-wide sequencing, reduced representation, and RNA sequencing, and NGS based technologies have not only supplied a plethora of information and source of markers. A transcriptome profile of fiber during its early stage in a fuzz less/lint less mutant of cotton plant revealed considerable deferential gene expression using a high-throughput Solexa genome analyzer [42, 45]. High throughput, genome-wide transcriptome study revealed that genes involved in the fiber elongation pathway were downregulated whereas, defense responsive genes were increased in cotton during drought stress [31, 47]. Fiber length ranging from 29 to 34 mm is thought to be long while above 34 mm ranks as extra-long and 22–28 mm lies as medium to medium long [5, 8, 9]. It is general consideration that longer fibers are often thought to be finer and stronger in comparison with short length fiber which are prime requirement of spinning process i.e., opening, cleaning, ginning, carding, drafting and combing [15]. On contrary, short fiber length may cause increased twisting during spinning resulting poor quality yarn [3, 10, 22].
Long-chain fatty acids (VLCFAs) promote fiber development by stimulating ethylene production leading to pectin biosynthesis and scaffold building; Ca2+-dependent protein kinase, reactive oxygen species (ROS), and sucrose synthase also aid in fiber cell elongation However, hypotheses, however, do not explain the molecular mechanism behind natural genetic changes in fiber length between genotypes. Fibers are categorized according to their fineness or softness which is basically the measurement of fineness and maturity called micronaire [36]. Fine fiber yarn has more fibers per cross- section resulting stronger and higher-quality yarn, premium range of fiber micronaire is from 3.7–4.2 μg/inch [49]. Several research have been undertaken to identify candidate genes involved in fiber growth, such as fiber micronaire improvement. Bajwa et al. [4] confirmed that the GhEXPA8 fiber expansion gene plays a vital role in extending fiber length and improving fiber micronaire by transforming GhEXPA8 in a cotton variety NIAB-846. The fiber strength of cotton is critical therefore, fiber with 30 g/tex are considered as stronger fiber while 23.5–25.4 g/tex are base line fiber [7]. In addition, fiber strength depends upon the deposition of cellulose inside the fibe, longer fiber contains more cellulose resulting stronger and good quality yarn [15, 21].
Lu et al. [30] has described 2200 common genes expressed differentially in high quality fiber chromosomal substitution lines but not found in low quality fiber lines. These genes were found to be involved in a range of metabolic processes regulating fiber strength. The genes XLOC075372 (snakin-1), XLOC036333 (mannosyl-oligosaccharide-a-mannosidase (MNS1), and XLOC029945 (FLA8) were discovered controlling regulation of cotton fiber strength which may need further investigation. Cotton fiber cell wall synthesis and strength might be regulated by RLKs, suggesting it a potential genes facilitating collaboration of both SCW assembly and fiber strength in other plants [20, 30]. The purpose of this study was to examine potential genes associated with various stages of fiber development, such as initiation, elongation, and primary and secondary cell wall formation during 0, 5, 10, 15, 20 and 25 DPA between very short fiber wild variety Gossypium darwinii; darwinii 5–7 and extra-long fiber sea-island Gossypium barbadense; XH-18 through transcriptome sequencing analysis. Furthermore, discovering DEGs associated to fiber development and determining the impact of wild allotetraploid cotton species i.e., Gossypium darwinii and Gossypium barbadense, on fiber quality. The comparisons were carried out between RNA-seq data collected from different fiber growing phases of both darwinii5–7 and XH-18. GFOLD analysis revealed a large number of DEGs, and the accuracy of GO functional enrichment and KEGG metabolic pathway analysis was carried out to confirm the DEGs further validated by quantitative real-time PCR (qRT-PCR). The findings demonstrated that the candidate genes are strongly linked to fiber development making a significant contribution to the mechanism of cotton fiber initiation, elongation and SCW formation.
Material and methods
Plant material
Gossypium barbadense L., XH-18 (EG16002) is the highly stable cultivar with extra-long fiber characteristics. Gossypium darwinii G. Watt., darwinii 5–7 (EG16004) is an exotic material with good fiber strength, drought and salt tolerance, nectarilessness traits which can be maintained through strict self-pollination. The parental materials were developed for boll tagging and fiber sample collecting at China's National Wild Cotton Germplasm Nursery located at Hainan province.
Transcriptome analysis
RNA extraction, library preparation and RNA sequencing
The TRlzol reagent was used to isolate total RNA from all 36 collected samples according to the manufacturer's instructions (Life technologies, California, USA). The RNA integrity and concentration were determined using an Agilent 2100 Bioanalyzer (Agilent Technologies, Inc., Santa Clara, CA, USA). NEBNext Poly (A) mRNA Magnetic Isolation Module was used to extract the mRNA (NEB, E7490). The NEBNext Ultra RNA Library Prep Kit for Illumina (NEB, E7530) and NEBNext Multiplex Oligos for Illumina (NEB, E7530) were used to construct cDNA library as per manufacturer's recommendations (NEB, E7500). The enriched mRNA was fragmented to generate around 250–300 bp RNA inserts to prepare first and second-strand cDNA and the second cDNA. The End-repair/dA-tail and adaptor ligation were undertaken on double-stranded cDNA. Using Agincourt AMPure XP beads (Beckman Coulter, Inc.) and PCR, the necessary fragments were extracted and amplified. Finally, using an Illumina HiSeqX sequencing technology, the produced cDNA libraries of both Gossypium darwinii and Gossypium barbadense were sequenced.
Mapping reference genome-based readings for transcriptome analysis
Peal script was used to delete low quality reads such as just adapter, unknown nucleotides greater than 5%, or Q20 20% (percentage of sequences with sequencing error rates less than 1%). The Tophat2 software was used to map the clean reads that were filtered from the raw data to the genome of (Gossypium species) [23]. The aligners'aligned data in BAM/SAM format were further evaluated to rule out possible duplicate molecules. FPKM values were utilized to assess gene expression levels using the Cufflinks programmer (fragments per kilo base of exon per million fragments mapped) [43].
Differential gene expression analysis
DESeq and Q-value were used to estimate the difference in gene expression between Gossypium darwinii and Gossypium barbadense [2]. The FPKM value ratio was then used to calculate the differences in gene abundance between those samples. The false discovery rate (FDR) control technique was used to estimate the P-value threshold in a range of tests in order to calculate the significance of the differences. Only genes having an absolute value of log2 ratio ≥ 2 and an FDR significance score of < 0.01 were included in the analysis.
Annotation of sequences
BLASTX with cut-off E-value of 10–5 was used to match genes with numerous protein databases, including the non-redundant protein (Nr) database at National Center for Biotechnology Information (NCBI) and the Swiss-Prot database (SPD). Genes were also screened against the NCBI non-redundant nucleotide sequence (Nt) database using BLASTn at cut-off E-value of 10–5. The best BLAST hit (highest score) was adopted to search genes, as well as their protein functional annotation. The Nr BLAST findings were given to Blast2GO program for gene annotation by using gene ontology employing GO keywords [11]. This study mapped all annotated genes and counted the number of genes associated with each phrase. A Perl script was used to plot GO functional classification for the unigenes having a GO term to study the distribution of gene functions. TopGo (an R package) was used to enrich and refine the final annotation [1]. The gene sequences were also matched to the Clusters of Orthologous Groups (COG) database to predict and categorize their activities [41]. The assembled sequences were assigned to KEGG pathways by using Perl script.
Comparison of expression patterns of DEGs
The short time-series expression miner software (Carnegie Mellon University, USA) was used to determine the temporal expression patterns of genes controlling growth phases of cotton fiber. GO enrichment analysis was also performed to assess the expression patterns and gene function identification through STEM analysis [12].
Development of primers for validation of DEG’s analysis
We selected the common DEG’s of Gossypium darwinii and Gossypium barbadense for transcriptome profiling at different stages of fiber development stages. The preparation of RNA samples was done as mentioned in “RNA extraction, library construction and transcriptome sequencing” section. The synthesis of cDNA was undertaken by using aTranscript All-in-One First-Strand cDNA Synthesis SuperMix (TransStart Top Green qPCR SuperMix kit (Transgen Biotech, Beijing,China), on the ABI 7500 fast Real-Time PCR System (Applied Biosystems, USA). The qPCR was done via adopting protocol of TransStart Top Green qPCR SuperMix kit (Transgen Biotech, Beijing,China), on the ABI 7500 fast Real-Time PCR System (Applied Biosystems, USA).
The housekeeping gene i.e., β-Actin gene was kept as the reference to normalize the relative expression levels, with its primer sequences: F: 5′-ATCCTCCGTCTTGACCTTG-3, and R: 5′-TGTCCGTCAGGCAACTCAT-3′. Following conditions for qRT-PCR was kept; 20 μL systems, one cycle of 94 °C for 30 s; 40 cycles of 94 °C for 5 s, 60 °C for 34 s, and one cycle of 60 °C for 60 s. To validate the result of qRT-PCR, three biological and technical replicates were used. 2−ΔΔCt method was adopted for the quantification of relative gene expression level.
Results
Total 12 RNA-seq libraries were generated and sequenced individually. The low-quality reads were eliminated and approximately 818 million clean reads data were attained and average reads per sample were calculated 45.49 million. In an alignment analysis, the clean data mapped on genome of Gossypium darwini were 86.93–88.5% which was not less than 93.883% of the Q30 with GC content of 43.917%, showing the reliability and quality of RNA sequencing data. Likewise in Gossypium barbadense, after elimination of low-quality reads, the clean reads were found approximately 844 million with average reads per sample were 46.925 million.
Differences in phenotypes of fiber traits between G. darwinii and G. barbadense and in F1 Population
The fiber related traits i.e., fiber length (mm), fiber elongation (%) and fiber strength (cN/tex) of Gossypium darwini were found non-significant compared to Gossypium barbadense whereas, fiber traits i.e., micronaire (%) value was significantly higher in Gossypium darwini compared to Gossypium barbadense. These results indicated that both Gossypium darwinii and Gossypium barbadense were more closely related species (Table 1). Similar results were obtained in number of genes involved in development of fiber at 0, 5, 10, 15, 20 and 25 DPA showing less in number for Gossypium darwini and greater in number for Gossypium barbadense from the clean reads data (Tables 2 and 3). The phenotypic observations between Gossypium darwinii and Gossypium barbadense at final stage comapred to F1 generation propose that there were substantial differences in final fiber quality of both species. The number of genes which found associated to initiation and elongation of fiber shown remarkable gradual difference for fiber quality.
Identification of DEGs involved in fiber development
A sum of 12,748 DEG genes were identified of which 6 out of 18 libraries were expressing identified genes. Pearson correlation Coefficient (PCC) analysis was implemented to examine the expression of genes derived from 18 identified libraries for testing the correlation among experimental samples (Figure S1 and Figure S2). There were more than 90% similarities in gene expression of both line at same stage of fiber formation with exception of gene found between darwinii 5–7 and XH-18 at 10 and 25 DPA. These gene expression similarities were greater at the 10 and 20DPA stages whereas, 25DPA stage demonstrated a significant shift in gene expression during fiber development. The gene expression at 5DPA indicated a substantial modification in gene expression which was consistent at different stages of fiber development.
The FPKM value was employed to compute the gene expression quantity for discovering the DEGs responsible for fiber development. The Generalized Fold Change (GFOLD) algorithm was exerted to discover DEGs between Darwinii 5–7 and XH-18 at the same fiber development stage and discovered 12,784 genes of which 6,742 were downregulated and 6,402 were found upregulated. The FPKM values of these DEGs in the two lines were clustered, resulting dynamic shift in DEGs as displayed in heat map (Fig. 1, Figure S3 and Figure S4). High numbers of DEGs, approximately 99,220 were found in genotype XH-18; whereas, 75,307 DEGs were found in Gossypium darwinii. Surprisingly, all DPAs revealed a significant degree of similarity in DEGs expression pattern for both Darwinii 5–7 and XH-18.
Functional Annotations of DEGs
On the basis of molecular functions, biological processes and cellular components 32,539 identified DEGs were categorized into 49 GO terms which further helped for understanding their role in fiber development and putative function (Fig. 2).
In the biological and cellular processes, the majority of differentially expressed genes (DEGs) were classified under metabolic processes, comprising 3,419 genes (10.51%) and 3,132 genes (9.63%), respectively. Only 31 genes were associated with growth processes. In the cellular component category, the largest subgroups were cell and cell part (2,294 genes, 7.05%), followed by membrane (2,184 genes, 6.71%).
The molecular function group was further comprised of subgroups i.e., binding subgroup comprised of 3,488 genes making 10.72% and catalytic activity group was comprised of 3,120 genes making 9.59%. Attributing to the expression quantity in the recent RNA sequencing data, 31 genes were found common throughout different DPAs stages of darwinii and 377 genes were common in Gossypium barbadense XH-18.
Hyper geometric model distribution to do GO enrichment and KEGG pathway analysis via adopting R software was undertaken to investigate the potential function of identified DEGs, respectively. The enrichment analysis was significant at FDR ≤ 0.05. The annotation and classification of 14, 584 upregulated genes into 48 GO terms were performed according to their molecular function, biological process and cellular component (Fig. 3).
In the biological processes, majority of (1,511 genes, 10.36%) up-regulated DEGs were involved in metabolic process resulting cell growth, related to cellular process, 1,357 genes which makes 9.30% were found responsible for cell physiology, growth and maintenance. Under the cellular component category, the most abundant genes were associated with cell and cell part (1,052 genes, 7.21%), followed by membrane (1,034 genes, 7.08%), and membrane part (958 genes, 6.57%). In the molecular function category, the major subgroups included binding (1,507 genes, 10.33%), catalytic activity (1,360 genes, 9.33%),"structural molecule activity (124 genes, 0.85%), and transporter activity (118 genes, 0.81%). A total of 17,955 downregulated genes were annotated and grouped into 48 GO terms across three categories: biological process, molecular function, and cellular component (Fig. 4).
Within the biological process category, the majority of downregulated differentially expressed genes (DEGs) were involved in metabolic processes (1,908 genes, 10.62%), cellular processes (1,775 genes, 9.89%), and single-organism processes (910 genes, 5.06%). In the biological process category, the majority of down-regulated DEGs were part of metabolic process (1,908 genes, 10.62 percent of the 17,955 down-regulated genes), cellular process (1,775 genes, 9.89%), and single-organism process (910 genes, 5.06%). In the cellular component and molecular function categories, the most represented subgroups were binding (1,981 genes, 11.03%), cellular process (1,775 genes, 9.89%), catalytic activity (1,760 genes, 9.80%), and cell/cell part (1,242 genes, 6.92%) (Figs. 5, 6 and 7). A gradual increase in novel DEG expression was observed across different fiber developmental stages. A total of 175 novel genes common to G. darwinii and G. barbadense were identified during fiber development (Fig. 8). Overall, an antagonistic pattern was seen between up- and down-regulated genes, with up-regulated genes showing an increasing trend and down-regulated genes a decreasing one.
Differential gene expression analysis
DEGs were also involved in various metabolic activities, including hydrolase activity particularly (hydrolyzing O-glycosyl compounds), the polysaccharide and carbohydrate catabolic process lipid metabolism, cellulase activity, metabolic process, catalytic activity, and general catalytic functions. They played integral role in the degradation of cellulose, lignin and pectin. Additionally, pathways related to phenylalanine metabolism i.e., flavonoid, phenylpropanoid, stilbenoid, diarylheptanoid and gingerol biosynthesis as well as ubiquinone, and other terpenoid-quinone biosynthesis and the degradation of aromatic compounds were found to significantly contribute to fiber production.
Differentially expressed plant hormones
Phytohormones i.e., brasinosteroids, auxin/IAA, cytokinin, ethylene, gibberellic acid and abscisic acid were have been found to play an important role in the fiber development. The differentially expressed genes associated with these hormones were catageorized into functional groups were categorized into functional categories i.e., cell response or signal transduction, redox homeostasis, protein, and energy or carbohydrate metabolism.
Differentially expressed transcription factors
Several transcription factor including WRKY, ranBP-type, C2H2 and C3HC4-type zinc finger-comprising protein, MYBP, JAZ, NAC, CKS, CNKX1E, MYO, bHLH, HAP, GAUT, BUB linked with DEGs, were found contribute to fiber initiation and elongation in the allotetraploid Gossypium barbadense and Gossypium darwinii species. Additionally, expensins, arabinogalactan proteins (AGPs), tubulin, plant receptor-like kinases (RLKs), leucine rich repeats, LRR-family protein and signal transduction coding for mitogen-activated protein kinase (MAPK) cascade involved in signal transduction have been implicated in cell elongation and SCW production.
Validation of novel differentially expressed genes (DEGs)
Subsequently, GO term analysis and KEGG annotation identified 20 novel expressed DEGs, including those involved in putative ribosomal RNA methyltransferase activity, which may enhance cell proliferation during cell cycle, and in translation elongation factor activity, ATP and GTP binding, translational elongation, nucleotide binding (XLOC_100202, XLOC_067150). The rest of the gene/s identified were associated with transporterand phosphodiesterase activity along with lipid metabolic process, phosphoric diester hydrolase activity (XLOC_071210, XLOC_074002) and ligase activity, zinc ion binding, metal ion bindings (XLOC_102238, XLOC_077416), were randomly selected and commonly found in Gossypium barbadense and Gossypium darwinii on the basis of venn diagram analysis for the confirmation of validity of RNA-seq data and β-Actin was kept as house-keeping gene (Table 4). The results showed that 18 genes (90 percent of 20 genes) had extremely congruent with the expression patterns of RNA-seq data, whereas the remaining two genes (10%) had only partially consistent (Table 4, Fig. 9, 10 and 11).
Discussion
Upland cotton (Gossypium hirsutum L.) is most widely used textile fiber source around the world. Cotton breeders have long sought to improve fiber quality by employing various genetics based approaches. In current study, the fiber development characteristics and quality attributes of fiber development and quality traits of darwinii 5–7 and XH-18 used as recurrent parent are described and presented in Tables 1, 2 and 3. The F1 genotypes exhibited superior fiber length and strength of the F1 genotypes compared to the parents. A sum of 32,539 DEGs were identified, expressed 6 out of 18 transcriptome libraries. Pearson Correlation Coefficient (PCC) analysis was emploed to evlauate the correlations between the experimental samples by examining the gene expression profile across these 18 libraries (Figs. 1 and 2). Comparative analysis revealed that gene expressions patterns between the two lines at the same fiber developmental stages were over 90% similar, with the exception of samples at 0, 5, 10, 20, and 25 DPA. These differences suggest a major transition in gene expression at the 10 DPA stage, indicating it as a critical point in fiber development.
The methyltransferase and ATPase binding activities associated with gene XLOC_100202 (Nop2) were implicated in cell elongation through stages of protein translation pathway during fiber development. These functions, including RNA binding and methylation, have been represented in previous studies [32]. Conversely, the gene XLOC_067150, which encodes a V-type proton ATPase subunit E-like protein is involved in ATPase binding activity. This gene plays a critical role for acidifying and regulating the pH of intracellular compartments of cell and specific particular cell types [44]. XLOC_071210, a probable aquaporin PIP2-8 gene, exhibits molecular function related to water channel activity and biological functions linked to the response to abscisic acid, which contributes to cell development and growth. The gene displayed both up and down regulated expression patterns during fiber development of Gossypium darwinii and Gossypium barbadense, respectively as confirmed by qRT-PCR results. RNA transcriptomics analysis showed consistently high expression at 0, 5, 10, 15, 20 and 25 DPA (Fig. 13 and 14). The PIP2, an intrinsic plasma membrane protein, is known to undergo novel post-translational regulation in plant cells [39]. Additionally, four uncharacterized genes XLOC_074002, XLOC_102238, XLOC_077416 and XLOC_091797 exhibited high expression levels in RNA sequence data across all fiber development stages (0, 5, 10, 15, 20 and 25DPA). These findings were further validated through gene expression analysis by using qRT-PCR for both species of cotton such as Gossypium darwinii and Gossypium barbadense. In total, 250 DEGs were identified that are involved in flavonoid pathways activities during fiber development in Gossypium barbadense. Previous studies have highlighted the importance of flavonoid pathway-related genes in fiber development, especially in both colored and wild cotton varieties [27]. Numerous flavonoid-associated genes i.e., CHS, LDOX, ANS, and F3H showed high expression levels during early fiber development in both Sea Island (extra-long staple) and upland cotton. These genes are believed to play a role in regulating plant development. Furthermore, some studies suggest that the flavonoid pathway may have contributed to cotton fiber evolution through negative selection during artificial domestication.
DEGs were also discovered to play roles in various biological processes and enzymatic activities, including hydrolase activity, hydrolyzing O-glycosyl compounds, polysaccharide catabolic process, glycosyl bonds clevage, carbohydrate metabolic, cellulase activity, and cellulose catabolic process. DEGs also participated in phenylalanine metabolism as well as in the biosynthesis of flavonoids, phenylpropanoids, stilbenoids, diarylheptanoids, gingerols, ubiquinone, and other terpenoid-quinones. Moreover, the degradation of aromatic compounds has been shown to play a significant role in the fiber production [16, 18, 19, 34].Certain genes related to fiber development, such as GhKCS13 (Qin et al. 2007), GbEXPA2 (Li et al. 2015), and GbEXPATR [29] have been identified in previous studies. One notable gene, Beta-glucuronosyltransferase (GlcAT14 A) has been found to function as a membrane associated component involved in plant-type secondary cell wall biogenesis. It exhibits xylosyltransferase activity during the glucuronoxylan biosynthetic process and is localized in the golgi apparatus, where it contributes to the formation of galactosylgalactosylxylosyl protein 3-beta-glucuronosyl. In Arabidopsis, beta-glucuronosyltransferase (AtGlcAT14 A) plays an important role in the biosynthesis of type II arabinogalactan (AG). This enzyme belongs to the glycosyltransferase family 14 (GT14) of the Carbohydrate Active Enzyme (CAZy) database and is essential for the regulation of plant growth and normal cellular functions. When transient expressed in Nicotiana benthamiana, the protein localized to the golgi apparatus [25]. supporting its role in xylan biosynthesis, a major component of the SCW in mature fiber cells.
Conclusion
This study focuses on comparing the transcriptome profiles of two cotton parents, Gossypium darwinii (Darwinii5-7) and Gossypium barbadense (XH-18), to explore the molecular mechanisms underlying fiber development. It aims to uncover the contributions of both species to fiber quality, offering a valuable foundation for improving fiber yield and quality in cultivated cotton. Notably, F1 hybrids exhibited enhanced biomass during early vegetative growth. A genome-wide transcriptome analysis across fiber development stages (0–25 days post-anthesis) revealed numerous DEGs in both parents. Many DEGs were overexpressed in fibers of both species, with G. barbadense showing higher expression levels, particularly during fiber initiation and elongation stages. These findings suggest that Gossypium barbadense may possess a more active regulatory network related to cell elongation and fiber growth. The observed differences in fiber traits likely result from significant changes in genes controlling cell expansion and associated metabolic processes. Overall, the data indicate that early fiber cell development is largely influenced by coordinated gene expression involved in both cell growth and metabolism.
Data availability
The data can be made available on reasonable request. The data is deposited at NCBI Bio-project number PRJNA1159014 (https://www.ncbi.nlm.nih.gov/bioproject/?term = PRJNA1159014).
Abbreviations
- Xh:
-
Xinhai
- GD:
-
Gossypium darwinii
- GB:
-
Gossypium barbadense
- GC:
-
Guanine-Cytosine Content
- SCW:
-
Secondary Cell Wall
- DPA:
-
Days Post Anthesis
- FL:
-
Fiber Length
- FS:
-
Fiber Strength
- FM:
-
Fiber Micrinaire
- NGS:
-
Next Generation Sequencing
- BAM/SAM:
-
Binary Alignment Map/Sequence Alignment Map Format
- FPKM:
-
Fragments per kb per million of the mapped reads
- GFOLD:
-
Generalized Fold Change
- STEM:
-
Short Time-series Expression Miner
- GO:
-
Gene Ontology
- DEGs:
-
Differentially Expressed Genes
- KEGG:
-
Kyoto Encyclopedia of Genes and Genomes
- COG:
-
Clusters of Ortholog Groups
- FDR:
-
False Discovery Rate
References
Alexa A, Rahnenfuhrer J. topGO: Enrichment analysis for Gene Ontology. R package version. 2010;2:2010.
Anders S, Huber W. 2010. Differential expression analysis for sequence count data. Genome Biol. 11. https://doi.org/10.1186/gb-2010-11-10-r106
Azzouz B, Ben Hassen M, Sakli F. Adjustment of Cotton Fiber Length by the Statistical Normal Distribution: Application to Binary Blends. J Eng Fiber Fabr. 2008;3(3):155892500800300300. https://doi.org/10.1177/155892500800300304.
Bajwa KS, Shahid AA, Rao AQ, Bashir A, Aftab A, Husnain T. Stable transformation and expression of GhEXPA8 fiber expansin gene to improve fiber length and micronaire value in cotton. Front Plant Sci. 2015;6:838. https://doi.org/10.3389/fpls.2015.00838.
Bölek Y, Çokkızgın H, Bardak A. Genetic analysis of fiber traits in cotton. KSÜ Doğa Bil Derg. 2014;17(1):15–20.
Bourland FM, Gbur EE. Relationships of plant trichomes to yield and fiber quality parameters in upland cotton. J Cotton Sci. 2017;21(4):296–305.
Bradow JM, Davidonis GH. Quantitation of fiber quality and the cotton production-processing interface: A physiologist’s perspective. J Cotton Sci. 2000;4(1):34–64.
Chandnani R, Kim C, Guo H, Shehzad T, Wallace JG, He D, Zhang Z, Patel JD, Adhikari J, Khanal S, Paterson AH. Genetic Analysis of Gossypium Fiber Quality Traits in Reciprocal Advanced Backcross Populations. The Plant Genome. 2018;11(1):170057. https://doi.org/10.3835/plantgenome2017.06.0057.
Chen X, Guo W, Liu B, Zhang Y, Song X, Cheng Y, Zhang L, Zhang T. Molecular mechanisms of fiber differential development between G. barbadense and G. hirsutum revealed by genetical genomics. PLoS One. 2012;7(1):e30056. https://doi.org/10.1371/journal.pone.0030056.
Cole RA, Fowler JE. Polarized growth: maintaining focus on the tip. Curr Opin Plant Biol. 2006. https://doi.org/10.1016/j.pbi.2006.09.014.
Conesa A, Götz S, García-Gómez JM, Terol J, Talón M, Robles M. Blast2GO: A universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21(18):3674–6. https://doi.org/10.1093/bioinformatics/bti610.
Ernst J, Bar-Joseph Z. STEM: A tool for the analysis of short time series gene expression data. BMC Bioinformatics. 2006;7:1–1. https://doi.org/10.1186/1471-2105-7-191.
Fang L, Gong H, Hu Y, Liu C, Zhou B, Huang T, Wang Y, Chen S, Fang DD, Du X, Chen H, Chen J, Wang S, Wang Q, Wan Q, Liu B, Pan M, Chang L, Wu H, Mei G, Xiang D, Li X, Cai C, Zhu X, Chen ZJ, Han B, Chen X, Guo W, Zhang T, Huang X. Genomic insights into divergence and dual domestication of cultivated allotetraploid cottons. Genome Biology. 2017;18:1–13. https://doi.org/10.1186/s13059-017-1167-5.
Fang L, Tian R, Chen J, Wang S, Li X, Wang P, Zhang T. Transcriptomic analysis of fiber strength in upland cotton chromosome introgression lines carrying different Gossypium barbadense chromosomal segments. PLoS One. 2014;9(4):e94642. https://doi.org/10.1371/journal.pone.0094642.
Foulk J, Meredith W, McAlister D, Luke D. 2009. Fiber and yarn properties improve with new cotton cultivar. J. Cotton Sci. 13.
Gou JY, Wang LJ, Chen SP, Hu WL, Chen XY. Gene expression and metabolite profiles of cotton fiber during cell elongation and secondary cell wall synthesis. Cell research. 2007;17(5):422–34. https://doi.org/10.1038/sj.cr.7310150.
Guan X, Song Q, Chen ZJ. Polyploidy and small RNA regulation of cotton fiber development. Trends Plant Sci. 2014. https://doi.org/10.1016/j.tplants.2014.04.007.
Han LB, Li YB, Wang HY, Wu XM, Li CL, Luo M, Wu SJ, Kong ZS, Pei Y, Jiao GL, Xia GX. The dual functions of WLIM1a in cell elongation and secondary wall formation in developing cotton fibers. Plant Cell. 2013;25(11):4421–38. https://doi.org/10.1105/tpc.113.116970.
Hua S, Wang X, Yuan S, Shao M, Zhao X, Zhu S, Jiang L. Characterization of pigmentation and cellulose synthesis in colored cotton fibers. Crop Sci. 2007;47(4):1540–6. https://doi.org/10.2135/cropsci2006.12.0835.
Huhtanen P, Rinne M, Nousiainen J. Evaluation of the factors affecting silage intake of dairy cows: A revision of the relative silage dry-matter intake index. Animal. 2007;1(5):758–70. https://doi.org/10.1017/S175173110773673X.
Jamshed, M., Jia, F., Gong, J., Palanga, K.K., Shi, Y., Li, J., Shang, H., Liu, A., Chen, T., Zhang, Z., Cai, J., Ge, Q., Liu, Z., Lu, Q., Deng, X., Tan, Y., Rashid, H. or, Sarfraz, Z., Hassan, M., Gong, W., Yuan, Y., 2016. Identification of stable quantitative trait loci (QTLs) for fiber quality traits across multiple environments in Gossypium hirsutum recombinant inbred line population. BMC Genomics 17. https://doi.org/10.1186/s12864-016-2560-2
Ji SJ, Lu YC, Feng JX, Wei G, Li J, Shi Yong-Hui YH, Fu Q, Liu D, Luo JC, Zhu YX. Isolation and analyses of genes preferentially expressed during early cotton fiber development by subtractive PCR and cDNA array. Nucleic Acids Res. 2003;31(10):2534–43. https://doi.org/10.1093/nar/gkg358.
Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: Accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14:1–3. https://doi.org/10.1186/gb-2013-14-4-r36.
Kim HJ, Triplett BA. Cotton fiber growth in planta and in vitro. Models for plant cell elongation and cell wall biogenesis. Plant Physiol. 2001;127(4):1361–6. https://doi.org/10.1104/pp.010724.
Knoch, E., Dilokpimol, A., Tryfona, T., Poulsen, C.P., Xiong, G., Harholt, J., Petersen, B.L., Ulvskov, P., Hadi, M.Z., Kotake, T., Tsumuraya, Y., Pauly, M., Dupree, P., Geshi, N., 2013. A β-glucuronosyltransferase from Arabidopsis thaliana involved in biosynthesis of type II arabinogalactan has a role in cell elongation during seedling growth. Plant J. 76. https://doi.org/10.1111/tpj.12353
Li H, Dou L, Li W, Wang P, Zhao Q, Xi R, Pei X, Liu Y, Ren Z. Genome-wide identification and expression analysis of the Dof transcription factor gene family in Gossypium hirsutum L. Agronomy. 2018;8(9):186. https://doi.org/10.3390/agronomy8090186.
Li, P. tao, Wang, M., Lu, Q. wei, Ge, Q., Rashid, M.H. or, Liu, A. ying, Gong, J. wu, Shang, H. hong, Gong, W. kui, Li, J. wen, Song, W. wu, Guo, L. xue, Su, W., Li, S. qi, Guo, X. ping, Shi, Y. zhen, Yuan, Y. lu, 2017. Comparative transcriptome analysis of cotton fiber development of upland cotton (Gossypium hirsutum) and chromosome segment substitution lines from G. hirsutum × G. barbadense. BMC Genomics 18. https://doi.org/10.1186/s12864-017-4077-8
Li X, Yuan D, Zhang J, Lin Z, Zhang X. Genetic Mapping and Characteristics of Genes Specifically or Preferentially Expressed during Fiber Development in Cotton. PLoS One. 2013;8(1):e54444. https://doi.org/10.1371/journal.pone.0054444.
Li Y, Tu L, Pettolino FA, Ji S, Hao J, Yuan D, Deng F, Tan J, Hu H, Wang Q, Llewellyn DJ, Zhang X. GbEXPATR, a species-specific expansin, enhances cotton fibre elongation through cell wall restructuring. Plant Biotechnol J. 2016;14:951–63. https://doi.org/10.1111/pbi.12450.
Lu, Q., Shi, Y., Xiao, X., Li, P., Gong, J., Gong, W., Liu, A., Shang, H., Li, J., Ge, Q., Song, W., Li, S., Zhang, Z., Rashid, M.H., Peng, R., Yuan, Y., Huang, J., 2017. Transcriptome analysis suggests that chromosome introgression fragments from sea island cotton (Gossypium barbadense) increase fiber strength in upland cotton (Gossypium hirsutum). G3 Genes, Genomes, Genet. 7. https://doi.org/10.1534/g3.117.300108
Lu X, Chen X, Mu M, Wang J, Wang X, Wang D, Yin Z, Fan W, Wang S, Guo L, Ye W. Genome-wide analysis of long noncoding rnas and their responses to drought stress in cotton (gossypium hirsutum l.). PLoS One. 2016;11(6):e0156723. https://doi.org/10.1371/journal.pone.0156723.
Maekawa S, Yanagisawa S. Ribosome biogenesis factor oli2 and its interactor brx1–2 are associated with morphogenesis and lifespan extension in arabidopsis thaliana. Plant Biotechnol. 2021;38(1):117–25. https://doi.org/10.5511/plantbiotechnology.20.1224a.
Meinert MC, Delmer DP. Changes in Biochemical Composition of the Cell Wall of the Cotton Fiber During Development. Plant Physiol. 1977;59(6):1088–97. https://doi.org/10.1104/pp.59.6.1088.
Pang, C.Y., Wang, H., Pang, Y., Xu, C., Jiao, Y., Qin, Y.M., Western, T.L., Yu, S.X., Zhu, Y.X., 2010. Comparative proteomics indicates that biosynthesis of pectic precursors is important for cotton fiber and Arabidopsis root hair elongation. Mol. Cell. Proteomics 9. https://doi.org/10.1074/mcp.M110.000349
Qin H, Chen M, Yi X, Bie S, Zhang C, Zhang Y, Lan J, Meng Y, Yuan Y, Jiao C. Identification of associated SSR markers for yield component and fiber quality traits based on frame map and upland cotton collections. PLoS One. 2015;10(1):e0118073. https://doi.org/10.1371/journal.pone.0118073.
Said JI, Lin Z, Zhang X, Song M, Zhang J. A comprehensive meta QTL analysis for fiber quality, yield, yield related and morphological traits, drought tolerance, and disease resistance in tetraploid cotton. BMC Genomics. 2013;14:1–22. https://doi.org/10.1186/1471-2164-14-776.
Salih H, Gong W, He S, Mustafa NS, Du X. Comparative transcriptome analysis of TUCPs in Gossypium hirsutum Ligon-lintless-1 mutant and their proposed functions in cotton fiber development. Mol Genet Genomics. 2019;294:23–34. https://doi.org/10.1007/s00438-018-1482-x.
Salih H, Leng X, He SP, Jia YH, Gong WF, Du XM. Characterization of the early fiber development gene, Ligon-lintless 1 (Li1), using microarray. Plant Gene. 2016;6:59–66. https://doi.org/10.1016/j.plgene.2016.03.006.
Santoni V, Verdoucq L, Sommerer N, Vinh J, Pflieger D, Maurel C. Methylation of aquaporins in plant plasma membrane. Biochem J. 2006;400(1):189–97. https://doi.org/10.1042/BJ20060569.
Seagull, R.W., Oliveri, V., Murphy, K., Binder, A., Kothari, S., 2000. Cotton fiber growth and development 2. Changes in cell diameter and wall birefringence. J. Cotton Sci. 4.
Tatusov RL, Galperin MY, Natale DA, Koonin EV. The COG database: A tool for genome-scale analysis of protein functions and evolution. Nucleic Acids Res. 2000. https://doi.org/10.1093/nar/28.1.33.
Tiwari SC, Wilkins TA. Cotton (Gossypium hirsutum) seed trichomes expand via diffuse growing mechanism. Can J Bot. 1995;73(5):746–57. https://doi.org/10.1139/b95-081.
Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, Van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–5. https://doi.org/10.1038/nbt.1621.
Vasanthakumar T, Rubinstein JL. Structure and Roles of V-type ATPases. Trends Biochem Sci. 2020. https://doi.org/10.1016/j.tibs.2019.12.007.
Wang F, Xu Z, Sun R, Gong Y, Liu G, Zhang J, Wang L, Zhang C, Fan S, Zhang J. Genetic dissection of the introgressive genomic components from Gossypium barbadense L. that contribute to improved fiber quality in Gossypium hirsutum L. Mol Breed. 2013;32:547–62. https://doi.org/10.1007/s11032-013-9888-y.
Wang S, Wang JW, Yu N, Li CH, Luo B, Gou JY, Wang LJ, Chen XY. Control of plant trichome development by a cotton fiber MYB gene. Plant Cell. 2004;16(9):2323–34. https://doi.org/10.1105/tpc.104.024844.
Wang XC, Li Q, Jin X, Xiao GH, Liu GJ, Liu NJ, Qin YM. Quantitative proteomics and transcriptomics reveal key metabolic processes associated with cotton fiber initiation. J Proteomics. 2015;114:16–27. https://doi.org/10.1016/j.jprot.2014.10.022.
Xu Z, Kohel RJ, Song G, Cho J, Alabady M, Yu J, Koo P, Chu J, Yu S, Wilkins TA, Zhu Y, Yu JZ. Gene-rich islands for fiber development in the cotton genome. Genomics. 2008;92(3):173–83. https://doi.org/10.1016/j.ygeno.2008.05.010.
Zumba J, Rodgers J, Indest M, 2017. Fiber micronaire, fineness, and maturity predictions using NIR spectroscopy instruments on seed cotton and cotton fiber, in and outside the laboratory. J Cotton Sci. 21.
Acknowledgements
We are indebted to give appreciation to State Key Laboratory of Cotton Biology, Institute of Cotton Research (ICR), Chinese Academy of Agricultural Sciences (CAAS); Anyang, China for providing umbrella for research program. We express profound sense of reverence to the entire research team, friends, and any other person who ever contributed; we have deep gratitude for you so much.
Funding
This research program was financially sponsored by the National Natural Science Foundation of China (32171994, 32072023) and the National Key R&D Program of China ((2021YFE0101200), Pakistan Science Foundation (PSF/CRP/18 thProtocol (07)), Public Sector Development Programme under Ministry of Planning, Development & Special Initiatives (PSDP Project No. 829), Pakistan Science Foundation and National Science Linkages Program (PSF/NSLP-9234) and International Foundation For Science, Sweden and COMSTECH, Islamabad (IFS-C-6500–1).
Author information
Authors and Affiliations
Contributions
K.W., A.D., and F.L. designed the experiments. A.D. and X.C. conceived the experiments and analyzed the results. X.C. Collected the material and A.D. and M.S. carried out the lab experiments and contributed equally. A.D. carried out all computational analyses. F.L., M.S., Y.X., Y.H., Y.H., Z.Z. and A.I. participated in the mapping experiments. A.D., Z.A. and M.S. drafted the manuscript. B.W., B.Z. and M.K.R.K., S.F., N.B.S.A.S., A.M.R provided critical revision, literature review, editing and proof-reading of the article. K.W. and F.L. finalized the manuscript. All authors read and approved the final manuscript.
Corresponding authors
Ethics declarations
Ethics approval and consent to participate
Not applicable.
Consent for publication
Not applicable.
Competing interests
The authors declare no competing interests.
Additional information
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary Information
12870_2025_6697_MOESM1_ESM.docx
Supplementary Material 1. Figure S1. Pearson Correlation Coefficient analysis within darwinii 5-7. Figure S2. Pearson Correlation Coefficient analysis within XH-18.Figure S3. A heat map representing novel DEGs from 12 samples in three replicates in darwinii 5-7 at different stages of DPAs like 0DPA, 5DPA, 10DPA, 15DPA, 20DPA and 25DPA, respectively. Figure S4. A heat map of novel DEGs from 12 samples in three replicates in XH-18 at different stages of DPAs like 0DPA, 5DPA, 10DPA, 15DPA, 20DPA and 25DPA, respectively
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.
About this article
Cite this article
Ditta, A., Cai, X., Shehzad, M. et al. Transcriptome dynamics along with expression analysis of key genes involved in fiber development between Gossypium barbadense and Gossypium darwinii. BMC Plant Biol 25, 691 (2025). https://doi.org/10.1186/s12870-025-06697-2
Received:
Accepted:
Published:
DOI: https://doi.org/10.1186/s12870-025-06697-2