- Research article
- Open Access
LCM and RNA-seq analyses revealed roles of cell cycle and translational regulation and homoeolog expression bias in cotton fiber cell initiation
BMC Genomics volume 22, Article number: 309 (2021)
Cotton fibers provide a powerful model for studying cell differentiation and elongation. Each cotton fiber is a singular and elongated cell derived from epidermal-layer cells of a cotton seed. Efforts to understand this dramatic developmental shift have been impeded by the difficulty of separation between fiber and epidermal cells.
Here we employed laser-capture microdissection (LCM) to separate these cell types. RNA-seq analysis revealed transitional differences between fiber and epidermal-layer cells at 0 or 2 days post anthesis. Specifically, down-regulation of putative cell cycle genes was coupled with upregulation of ribosome biosynthesis and translation-related genes, which may suggest their respective roles in fiber cell initiation. Indeed, the amount of fibers in cultured ovules was increased by cell cycle progression inhibitor, Roscovitine, and decreased by ribosome biosynthesis inhibitor, Rbin-1. Moreover, subfunctionalization of homoeologs was pervasive in fiber and epidermal cells, with expression bias towards 10% more D than A homoeologs of cell cycle related genes and 40–50% more D than A homoeologs of ribosomal protein subunit genes. Key cell cycle regulators were predicted to be epialleles in allotetraploid cotton. MYB-transcription factor genes displayed expression divergence between fibers and ovules. Notably, many phytohormone-related genes were upregulated in ovules and down-regulated in fibers, suggesting spatial-temporal effects on fiber cell development.
Fiber cell initiation is accompanied by cell cycle arrest coupled with active ribosome biosynthesis, spatial-temporal regulation of phytohormones and MYB transcription factors, and homoeolog expression bias of cell cycle and ribosome biosynthesis genes. These valuable genomic resources and molecular insights will help develop breeding and biotechnological tools to improve cotton fiber production.
Cotton is the largest renewable source of textile fiber in the world and is an important oil crop [1, 2]. Most widely cultivated cottons are allotetraploid, consisting of two sets of chromosomes from different origins across the old and new world, with the “mother” (A-genome like species) from Africa and the “father” (D-genome like species) from the Americas . The polyploidization took place in the new world ~ 1.5 million years ago , and resulted in five extant species that diversified over 300,000–600,000 years [4, 5]. Over the last 8000 years, two of them, Upland (Gossypium hirsutum L.) and Pima (G. barbadense L.) cottons were independently domesticated in NW S. America and the Yucatan Peninsula of Mexico, respectively, under strong human selection, leading to the modern annualized crops . To date, Upland cotton dominates in ~ 94% of world cotton production, with the remainder (~ 5%) being Pima cotton for its superior fiber quality , and other diploid cottons adapted to certain growth environments .
Cotton fiber is seed hair, and each individual cotton fiber initiates from a single epidermal-like cell from the ovular epidermis, which elongates dramatically to reach up to ~ 6 cm in length, one of the longest singular cells in the plant kingdom . Each cotton seed has over 20,000 semi-synchronically developed fiber cells, representing 20–30% of cells in the epidermal layer of the seed [8, 9]. There are two types of fiber, lint and fuzz fibers [10, 11], and only the lint fibers are ginned and used in textiles. Lint fiber cells are initiated before or on the day of anthesis (0 Day Post Anthesis, DPA), which induces production and transportation of the phytohormones auxin and gibberellin to promote fiber cell elongation . Fuzz fiber is produced approximately 5 days after the development of lint fibers and is much shorter than the lint fiber. Fiber cells undergo overlapping stages of elongation, cellulose biosynthesis, and maturation . During elongation phase, primary cell wall is synthesized, and the cell extends through some combination of anisotropic expansion and tip growth . Once fibers have reached full length, they transition to the stage of cellulose and secondary cell wall biosynthesis, which is accompanied by the maturation phase where fibers, now comprising ~ 95% pure cellulose, become metabolically inactive, leading to cell death, which happens before the boll opens .
In developing fiber cells, cell cycle is arrested at G1 or S phase . Cell differentiation and cell cycle regulation are linked in plants , which involves cyclin-dependent kinases (CDKs) , phytohormones (auxin and gibberellins), and key transcription factors. Cell cycle exit without mitotic division often leads to endoreduplication, which is common in eudicot plants . Previous studies showed that DNA content was increased by endoreplication during early fiber development [18, 19], which is accompanied with small RNA and DNA methylation changes in early fiber cells [20,21,22]. In human cancer cells, tumor suppressor p53 has been shown to induce cell-cycle arrest in response to impaired ribosome biogenesis , and vice versa, cell cycle arrest during tumor metastasis is fueled by upregulation of ribosome biogenesis .
Phytohormones and transcription factors also play crucial roles in fiber cell development [2, 25,26,27]. Auxin, brassinosteroids, ethylene, and gibberellins promote or induce fiber cell development [28,29,30], whereas cytokinin and abscisic acid inhibit fiber cell development . Previous studies have shown that auxin accumulates in ovules before fertilization, peaking around 2–3 DPA, with the level gradually decreasing through the elongation stage , while abscisic acid gradually accumulates in fiber cells during the initiation and elongation stages, with a peak around 10 DPA, followed by a decrease in later stages .
MYB (myeloblastosis) transcription factors are required for development of leaf trichomes and seed hair or fiber cells [33,34,35,36]. Several MYB transcription factor families are shown to regulate fiber cell initiation and elongation [34,35,36,37], while some MYB transcription factor members do not have an obvious role in cotton fiber cell development .
Regulatory networks involving phytohormones and transcription factors during early stages of fiber cell development remain poorly understood, partly because it is technically difficult to isolate fiber and epidermal cells during fiber cell initiation. Here we performed RNA-seq analysis using the fiber and epidermal cells that were separated by Laser Capture Microdissection (LCM) in early stages (0 and 2 DPA). We compared transcriptomes of LCM samples between the wild-type and fiber-less mutant (both in Upland background) and between Upland and Pima cottons, which revealed expression differences in the genes associated with cell cycle regulation and ribosome biosynthesis. Using pharmacological inhibitors for these pathways, we showed that cell cycle regulation and ribosome biogenesis affect fiber cell initiation and elongation in the ovules cultured in vitro. Moreover, phytohormone related genes were upregulated in the ovules and down-regulated in early fibers, suggesting a spatial-temporal role in fiber cell development. MYB-transcription factor related genes displayed expression divergence between fibers and ovules, while key cell cycle regulators were predicted to be epialleles . Together, these results support a new role for cell cycle and translational regulation in fiber cell initiation and elongation, which is accompanied by spatial-temporal regulation of transcriptional and phytohormonal networks.
Fiber cell initiation in GhTM-1, Gb3–79, and GhMD17
Our study involved two cultivated cottons, G. hirsutum L. acc. TM-1 (GhTM-1) and G. barbadense L. acc. 3–79 (Gb3–79), and the fiber-less mutant GhMD17 (G. hirsutum background). MD17 is a naturally occurring mutant with naked seed (lintless and fuzz fiberless) phenotype, which is proposed to be controlled by three genetic loci N1N1n2n2n3n3 . N1 locus is controlled by MYB25-like (aka, MML1_A12) gene on A12 chromosome, which generates small interfering RNAs from its 3′ end that post-transcriptionally silence both homoeologs of this gene , and n2 locus is related to MYB25-like gene on homoeologous D12 chromosome [37, 42]. These three cotton lines were chosen on the basis of their distinct fiber characteristics: GhTM-1 has both lint and fuzz fibers, Gb3–79 has almost exclusively lint fibers, and GhMD17 is fiberless (Fig. 1a). Both GhTM-1 and Gb3–79 have similar lint fiber properties in mature seeds, except for the lack of fuzz fiber on the Gb3–79 seed . The fuzz fiber varies by genotypes and present in the Upland cotton (G. hirsutum) but largely absent in Pima cotton (G. barbadense). The genetic basis of this difference is not fully understood and is likely related to mutation of MYB25-like on homoeologous D12 chromosome .
During fiber initiation stage (0 and 2 DPA), both GhTM-1 and Gb3–79 had fibers of ~ 10 μm in length at 0 DPA and ~ 50 μm at 2 DPA, with little difference between GhTM-1 and Gb3–79, although fiber growth was slightly faster in Gb3–79 than in GhTM-1 (Fig. 1c). GhMD17 lacked any visible fiber cell initials during those stages (Fig. 1b). Thus, GhTM-1 and Gb3–79 have a similar fiber initiation process, which is blocked in GhMD17.
RNA-seq analysis of fiber and epidermal cells separated by LCM
In wild-type ovules at 0 DPA and fiber-less mutants at both 0 and 2 DPA, LCM was used to cut the top layer of epidermal cells. For wild-type ovules at 2 DPA, LCM was used to cut along the epidermal surface to collect protruding fibers, followed by a second cut to harvest underlaying epidermal cells (Fig. 1d and Supplemental Figure 1a). LCM materials in each sample with two biological replicates were used for total RNA extraction (Supplemental Figure 1b), library construction, and sequencing. Sequencing reads were mapped onto the G. hirsutism genome v2.0 (GhTM-1 and GhMD17) and G. barbadense genome v1.0 (Gb3–79)  (Supplemental Table 1), and read counts were normalized to fragments per kilobase per million (FPKM) for gene expression analysis (see Methods).
Principal component analysis (PCA) showed a clear separation of gene expression variance by tissue types (PC1: 37% of variance), followed by species (PC2: 23% of variance) (Fig. 2a). Epidermal samples between the 0 and 2 DPA stages were not separated. We also determined the Pearson’s correlation coefficient of all 8 samples (averaged values of two biological replicates) and performed hierarchical clustering analysis among these samples (Supplemental Figure 2a). This revealed three major groups: fiber samples (GhTM-1 and Gb3–79 at 2 DPA), mutant epidermal samples (GhMD17 at 0 and 2 DPA), and Gb3–79 epidermal samples (0 and 2 DPA). GhTM-1 epidermal sample at 0 DPA was grouped with Gb3–79 epidermal samples, while GhTM-1 epidermal sample at 2 DPA and mutant epidermal samples were in the same group. These results indicate that gene expression differences between tissue types and between species of LCM samples were suitable for this study.
At 0 DPA, fiber cell initials develop in the wild type but not in the mutant. To identify differentially expressed genes in fibers at 0 DPA, we compared expression differences in the epidermal samples between the wild type (GhTM-1 and Gb3–79) and mutant (GhMD17). We identified 1295 and 1500 upregulated genes in GhTM-1 and Gb3–79 fiber cell initials, respectively, 981 of which (> 60%) were shared in both, while 314 and 519 were specific to GhTM-1 and Gb3–79, respectively (Fig. 2b, Supplemental Tables 2 and 3). In addition, 2101 and 2582 genes were down-regulated in GhTM-1 and Gb3–79 fiber cell initials, respectively, 1528 of which (> 60%) were shared in both (Supplemental Figure 2b, Supplemental Table 2).
In epidermal cells, there were 291 and 181 upregulated genes in GhTM-1 and Gb3–79, respectively, at 0 DPA (Supplemental Figure 3a and Supplemental Table 3a), and 583 and 409 upregulated genes in GhTM-1 and Gb3–79, respectively, at 2 DPA (Supplemental Figure 3b and Supplemental Table 3a). The number of shared genes between GhTM-1 and Gb3–79 was relatively small, only 14 genes (~ 5%) at 2 DPA and 121 genes (~ 20%) genes at 2 DPA, suggesting expression divergence in epidermal layer cells between Upland and Pima cotton. GO analysis showed overrepresentation of upregulated genes at 2 DPA in lipid metabolic process (p-value 0.00045) and DNA replication (p-value 2.7e-12), consistent with biological processes involved ovule and seed development.
At 2 DPA, we compared gene expression of LCM fiber cell with epidermal cell tissues within the same genotype (GhTM-1 or Gb3–79) and between the fiber cell (GhTM-1 and Gb3–79) and epidermal cell tissues in the fiber-less mutant (GhMD17). We identified 1861 and 2046 up-regulated genes in GhTM-1 and Gb3–79 fiber cells, respectively, 743 of which (> 35%) were shared in both, while 265 and 450 were unique to GhTM-1 and Gb3–79, respectively (Fig. 2c, Supplemental Tables 3 and 4). We also identified 2978 and 4269 down-regulated genes in GhTM-1 and Gb3–79 fiber cells, respectively, and 1626 of which (> 35%) were shared in both (Supplemental Figure 2c, Supplemental Table 4).
To exclude a possibility of the fiberless mutant (MD17) effect, we generated data between GhTM-1 and Gb3–79 in fiber and epidermal samples at 2DPA (Fig. 3 and Supplemental Table 4a). A total of 1259 and 1631 genes in fiber cells at 2 DPA were upregulated in GhTM-1 and Gb3–79, respectively (Fig. 3a), while 2764 and 3111 genes in epidermal cells were upregulated in GhTM-1 and Gb3–79, respectively (Fig. 3b). Among them, more upregulated genes were shared between GhTM-1 and Gb3–79 in epidermal layer (2131, > 68.5%) than in fiber (931, > 56.8%), which were consistent with the results including GhMD17 epidermal cells. These differentially expressed genes in the fiber cell tissues were used for further analysis.
Dramatic gene expression changes in early stages of cotton fiber cell development
Gene Ontology (GO) analysis of these differentially expressed genes showed a consistent overrepresentation of translation (ribosome activity), gene expression, and peptide, lipid, and fatty acid metabolic processes in the upregulated genes in fibers (Fig. 4 and Supplemental Table 4b). These results may indicate that ribosome biosynthesis, required for protein synthesis as well as lipid and fatty acid metabolism, were active during fiber cell initiation. Moreover, the upregulated genes between the species were overrepresented in ATP and purine metabolic processes in GhTM-1, and ADP and pyridine (similar structure of pyrimidine) metabolic processes were overrepresented in Gb3–79. ATP is a critical energy resource for fiber growth. In cultured cotton ovules, fibers release ATP as they grow, and when the ectoapyrase activity is blocked by the addition of polyclonal anti-apyrase antibodies or small molecule inhibitors, the medium ATP level rises and fiber growth is suppressed; low concentrations of hydrolyzable nucleotides ATPgammaS/ADPbetaS stimulate fiber growth . Thus, the difference of ATP and ADP metabolic processes between GhTM-1 and Gb3–79 may contribute to longer fibers in Gb3–79.
Among down-regulated genes in fibers, they were overrepresented in response to hormones and glucose metabolic process at 0 DPA and DNA replication and the cell cycle at 2 DPA. Together, these data may suggest that cell cycle and DNA replication have attenuated in fiber cells, which require active ribosome biosynthesis for translational regulation in fiber cells. Down-regulation of phytohormone-related genes in the fiber cells may suggest upregulation of these genes in ovules to support fiber elongation.
Pervasive homoeolog expression bias in fiber and epidermal cells
Allotetraploid cottons consist of A and D subgenomes (A and D) with 52,514 orthogroups representing 15,220 pairs of A and D homeologs in G. hirsutum and G. barbadense (Supplemental Tables 5 and 5a). Among them, ~ 20% of homeolog pairs were differentially expressed (Supplemental Figure 4a), with 10% more D homoeologs (1441-1550) expressed at higher levels (D > A) than the A homoeologs (1250-1402, A > D) in both fiber and epidermal cells (Supplemental Figure 4b,c). These data are consistent with the notion of dramatic gene expression divergence between the homoeologs, despite gene numbers and synteny are conserved between the A and D subgenomes in the allotetraploid cottons . This trend of biased expression is shared by 251 and 81 homoeolog pairs that encode cell cycle regulation genes and ribosomal protein subunit genes, respectively, with expression bias for ~ 20% of homeolog pairs involved in cell cycle regulation and for ~ 15% of homoeologs encoding ribosomal protein subunits (Supplemental Figure 4d). Interestingly, for cell cycle regulation genes, 10% more D than A homoeologs had biased expression, whereas this number of homoeologs with expression bias was increased to 40–50% more in D than in A homoeologs of ribosomal protein subunit genes (Supplemental Figure 4d). This high number of homoeologs with expression bias in the early stages of fiber cells may indicate that deep subfunctionalization of homoeologous genes during active developmental process of fiber cell initiation and elongation.
Translational regulation in early cotton fiber cell development
During fiber cell development, cell cycle is predicted to attenuate in 20–30% of epidermal cells, and those cells become fiber cell initials. Cell cycle arrest in a cell often leads to active ribosome biosynthesis . This is consistent with developmental process of fiber cell initiation, when the arrest of the cell cycle in the epidermal cells induces active ribosome biosynthesis for fiber cell initiation and elongation.
Arabidopsis has 231 annotated genes for cell cycle regulation, and cotton has a total of 845 putative cell-cycle regulation genes, including 366 involved in G1, 283 in G2, 122 in S, and 74 in M stages (Supplemental Table 6). Approximately 85% of all putative cell-cycle regulation genes were down-regulated in fiber cells (Fig. 5a and Supplemental Figure 5a). Cotton has 681 annotated genes encoding for ribosome protein subunits, consisting of large subunit, RPL and small subunit RPS, respectively (Supplemental Table 6). These multiple putative ribosomal subunit genes (618 genes) are grouped into 121 putative ribosomal proteins. In contrast to down-regulation of cell cycle related genes, approximately 75% of putative ribosomal protein genes were up-regulated, while 8% were down-regulated (Supplemental Figure 5b). Plants have more ribosomal protein genes than animals or yeasts , and many of those paralogues are differentially expressed during development and in UV-B responses . For example, in Arabidopsis RPL16A and RPL16B are expressed in developing pollen and early phase of lateral root initiation , a process similar to fiber cell initiation. AtRPS5A gene is strongly expressed in dividing cells, whereas AtRPS5B expression correlates with cell differentiation rather than cell division . These results suggest that during fiber cell initiation, cell cycle was attenuated, and ribosomal protein biosynthesis is activated and may regulate tissue-specific expression in fiber and epidermal cells.
Cell cycle is regulated by Target of Rapamycin Complex1 (TORC1), which consists of three proteins TARGET OF RAPAMYCIN (TOR), RAPTOR, and LETHAL WITH SEC THIRTEEN PROTEIN 8 (LST8) [49, 50]. ErbB-3 EPIDERMAL GROWTH FACTOR BINDING PROTEIN1 (EBP1) is a TOR target. TORC1 is also a key regulator of many other biological processes such as DNA and RNA synthesis and autophagy and is conserved among plants and animals . In cotton, TOR, RAPATOR, LST8, and EBP1 have multiple members, two homoeolog pairs each for TOR, RAPTOR, and EBP1, and one homoeolog pair for LST8. Expression levels among multiple members of each gene were similar and averaged for the analysis (Supplemental Table 7). In fiber cells at 2 DPA, TOR was down-regulated in both GhTM-1 and Gb3–79 (Fig. 5b), while RAPTOR and LST8 expression levels were similar to epidermal cells, and EBP1 was slightly up-regulated. In Arabidopsis, TOR is expressed in the primary meristem, embryo, and endosperm, but not in differentiated cells , while the mutant is embryonically lethal. Notably, Arabidopsis EBP1, a homolog of human EBP1 for ribosomal protein biosynthesis, is required for expression of cell cycle genes and regulates cell proliferation and size in plants . Thus, in fiber cells down-regulation of TOR may inhibit the cell cycle, and strong expression of EBP1 could promote ribosome biosynthesis for fiber cell expansion. Interestingly, both cell cycle regulators TOR and RAPTOR have been identified as epialleles, which are demethylated in the cultivated cotton species but not in the wild relatives .
To further test roles of ribosome biosynthesis and cell cycle regulation in fiber cell development, we used pharmacological assays in cultured GhTM-1 and GhMD17 ovules. Ribozinoindoles (Rbin-1) are potent and reversible triazinoindole-based inhibitors of eukaryotic ribosome biogenesis , and Roscovitine is a cyclin-dependent kinase inhibitor for cell cycles . Cotton ovules were isolated from fertilized flowers of GhTM-1 and GhMD17 at 0 DPA and incubated in BT medium for 21 days following a published protocol , with or without inhibitors. Inhibiting cell cycle by Roscovitine (2 μg/L) increased total fiber area and length at statistically significant levels, whereas impeding ribosome biosynthesis by Rbin-1 (40 μg/L) reduced the fiber area (Fig. 5c, d), and the ovule size was not affected by either treatment (Fig. 5d). Indeed, chemical inhibitors had no effects on ovule development in the GhMD17 that is fiberless (Fig. 5c left). The altered fiber area could result from increased fiber length, number, or both. To discern this, we further measured the number or density of fiber cell initials after 2 days of each treatment (Supplemental Figure 6). The result showed a similar trend of increasing and reducing fiber cell density by Roscovitine and Rbin-1 treatments, respectively. These results support the hypothesis that cell cycle attenuation and active ribosome biosynthesis promote fiber cell initiation and elongation.
Dynamic roles of phytohormones and MYB transcription factors in early fiber cells
Phytohormones play important roles in fiber cell development [2, 25], but the temporal regulation of gene expression in early fiber developmental stages is poorly understood. Using LCM, we unexpectedly uncovered genome-wide down-regulation of 43 phytohormone related genes with only a few up-regulated genes in early fiber cells (Fig. 6a and Supplemental Table 8). This is probably because most genes are upregulated in the epidermal layer cells that support fiber cells. Abscisic acid (ABA) is known to inhibit fiber cell initiation and elongation . ABA in ovules at 0 DPA accumulates at much higher levels in the Ligon-lintless 1 mutant than in the wild type . This is consistent with the down-regulation of key regulators in the ABA pathway (Fig. 6a), including phosphate type-2C HAB1 (a homeolog of ABA INSENSITIVE 1, ABI1; and ABI2)  and ABA receptor gene PYRABACTIN RESISTANCE-LIKE PROTEIN 9 (PYL9) (known to interact with ABI1 and ABI2) . Down-regulation of the ABA regulator and receptor genes may promote fiber cell initiation and elongation.
Gibberellins promote fiber cell elongation and increase ABA and auxin levels during fiber development to enhance fiber quality (strength) . Over-expression of gibberellin biosynthesis gene GIBBERELLIN 20-OXIDASE1 (GA20OX1) leads to more and longer fibers , while their homeologs GA20OX2 and GA20OX3 are largely expressed in the ovule. Consistently, GA20OX2 was upregulated in the ovule and down-regulated in fibers (Fig. 6a).
Ethylene regulates fruit ripening and inhibits vegetative growth [62, 63], and is reported to promote fiber elongation in cotton, which is mediated by very-long-chain fatty acid signaling pathway . However, the mechanism for ethylene in fiber initiation remains unclear. We found down-regulating of the ethylene signaling pathway genes like ETHYLENE INSENSITIVE 2 (EIN2) and EIN3 in fibers, which may suggest a negative role for ethylene in fiber cell initiation (Fig. 6a).
Auxin promotes fiber initiation and accumulates in ovules prior to fertilization . Auxin is transported from ovule to fiber cells during fiber cell development . Consistent with this notion, most genes encoding AUXIN RESPONSE FACTORs (ARFs) and PIN-FORMED6 (PIN6), a key auxin transporter, were upregulated in the ovules and down-regulated in early fiber cells (Fig. 6a). This is reminiscent of the finding that ARF2, ARF18 and PIN6 were expressed in ovules in the fiber initiation stage . Interestingly, although the majority of among ARF gene family members are upregulated in the ovules, a few such as ARF19a and ARF19b homoeolog pair, showed opposite expression patterns in the fiber cells (Supplemental Figure 7), suggesting functional divergence between homoeologs for fiber cell development in polyploid cotton.
Brassinosteroids promote both fiber cell initiation and elongation [28, 64]. A key regulator of brassinosteroid signaling, BRASSINAZOLE-RESISTANT 1 (BZR1), induces cotton fiber development . BSU1 is a dominant suppressor of BRI1 and is expressed in elongating cells in Arabidopsis . In cotton fiber cells, BSU1 was upregulated, whereas BZR1 and its semidominant suppressor gene (BES1)  were down-regulated (Fig. 6a), suggesting a role for brassinosteroids in fiber cell initiation.
In addition to phytohormones, MYB transcription factors regulate fiber cell development in seeds and trichome development in leaves [37, 38, 68]. Among 18 MYB transcription factor genes examined, nine were upregulated, and the other nine were down-regulated in early fiber cells (Fig. 6b and Supplemental Table 8). The upregulated genes included GhMYB25, GhMYB109, and cotton homologs of Arabidopsis MYB4, MYB66, MYB96, and MYB106 (Supplemental Table 8). Notably, cotton homologs of MYB106 (GhMYB25 in cotton) were up-regulated, while cotton MYB16 homologs were down-regulated in fibers (Fig. 6b), while Arabidopsis MYB16 and MYB106 are in the same clade (Supplemental Figure 8). Similarly, cotton MYB4 was upregulated, and MYB6 was down-regulated in early fiber cells (Fig. 6b and Supplemental Figure 8), but they are closed related MYB genes. These results may suggest functional divergence between closely related homologous genes during fiber cell development.
In this study, we used LCM to successfully dissect fiber and epidermal cells from the cotton ovules during fiber cell initiation and investigated gene expression changes between fiber cell initials and epidermal layer cells. Our data support a new model of fiber cell initiation and elongation, involving regulation of cell cycle and ribosome biosynthesis activities, as well as phytohormones and MYB transcription factors (Fig. 7).
Fiber cells are initiated in 20–30% of epidermal cells at or prior to 0 DPA and before fertilization. Fiber cell initials differentiate from epidermal layer cells and involve a rapid and dramatic transition from actively dividing ovular cells to elongation of new singular fiber cells. At the cellular and molecular levels, cell cycle is predicted to arrest and coupled with active biosynthesis of ribosome for translational activities including protein synthesis. This rapid process of cell differentiation is similar to the origin of human tumor cells, where cell cycle arrest and ribosome biosynthesis are intrinsically correlated [23, 24]. The tumor suppressor p53 has been shown to induce cell-cycle arrest in response to active ribosome biogenesis for rapid cell growth , and the transition associated with development and tumor metastasis is accelerated by upregulation of ribosome biogenesis during G1/S arrest . By this analogy, although the cause and consequence are unknown, some unknown factor(s) may trigger cell cycle arrest in pre-fiber cells to induce active ribosome biosynthesis to prepare for rapid cell initiation and elongation. Considering that each fiber cell does not divide, abundant ribosomes may promote translational regulation and protein synthesis, leading to rapid fiber cell elongation.
The TOR signaling pathway integrates nutrient, hormone, growth factor and environmental inputs to control cell cycle, proliferation, metabolism, and growth in diverse organisms, including plants . Loss of TOR function results in cell cycle arrest in G1 phase in human and yeast , which is consistent with the observation that cotton fibers cells are arrested at G1 stage [14, 70]. We found in early stages of fiber cells, TOR expression is downregulated (Fig. 5), which may induce ribosome biosynthesis and translation rate through the TOR signaling pathway . Moreover, EBP1 is upregulated in the cotton fiber cells, which can positively regulate ribosome biosynthesis as part of the TOR signaling pathway and auxin regulation in plants . It is likely that the TOR signaling pathway may regulate cell cycle and ribosome biosynthesis in the ovule epidermis for fiber cell initiation as well as elongation in early fiber cells.
Phytohormones are essential for fiber cell initiation and development [2, 25, 56]. While some phytohormones such as ABA and cytokinin have a negative role in fiber cell development, many others such as auxin, gibberellins, and brassinosteroids positively regulate fiber cell initiation and development. It is unexpected to find genome-wide down-regulation of phytohormone related genes such as ARFs, PIN6, BRZ1, EIN2, EIN3, and GA20OX2 in early fiber cells. This result may suggest an important role for spatial and temporal action of phytohormones in fiber cell development, which is revealed by the LCM approach in this study. Many phytohormone related genes such as ARFs and PIN6 are strongly expressed in epidermal layer cells  (Fig. 6). Auxin may actively accumulate in the ovules and is transported into young fiber cells, reflecting a spatial and temporal role of auxin in fiber cell initiation. Indeed, expression of cotton ARF2:GUS in Arabidopsis is localized in the trichome base but not in the body or tip , consistent with its spatial and temporal role in leaf trichome development as in cotton fiber cell initiation. Alternatively, there may be crosstalk between phytohormone pathways and other factors such as cell cycle regulation and ribosome biosynthesis to regulate fiber cell development. For examples, inhibiting TOR signaling pathway decreases ABA synthesis , and auxin activates TOR signaling pathway and promotes translational regulation  and some other associated pathways .
MYB transcription factors are required for the development of leaf trichomes as well as seed hair or fiber cells [33,34,35,36]. Cotton has two major groups of MYB transcription factor genes, one expressed in trichomes or leaves and another largely in ovules and fiber cells . Up and down-regulation of MYB transcription factor related genes may reflect different roles of these genes in trichome and fiber cell development in cotton. In Arabidopsis, MYB16 and MYB106 (GhMYB25 in cotton) are closely related genes that promote cuticle formation in the epidermal cells including trichomes  and promote both trichome branching and expansion , while the cotton homolog GhMyb25 regulates cotton fiber cell development . Thus, upregulation of MYB106 and down-regulation of MYB16 in early fiber cells may reflect their functional divergence in trichome and fiber cell development. This expression divergence between closely related homoeologs is common to many other genes including those in phytohormone, cell cycle, and ribosome biosynthesis pathways. The mechanism for this expression divergence may involve small RNAs including siRNAs and miRNAs, DNA methylation, and other epigenetic modifications [20,21,22, 34].
Finally, the genes involved in photoperiodic flowering and cell cycle regulation such as TOR and RAPTOR are reported to be epialleles, which are demethylated and highly expressed in the cultivated cotton species but not in the wild relatives . This result underscores important role for epigenetic regulation in selection and domestication of agronomic traits including fiber cell development, possibly fiber yield and quality. Our study has provided valuable genomic resources and new insight into roles of cell cycle regulation, ribosome biosynthesis, and spatial-temporal action of phytohormones and MYB transcription factors in early stages of fiber cell development. Understanding this complex regulation during fiber cell initiation and elongation may help us develop new breeding and biotechnological tools to improve cotton fiber yield and quality.
Two cultivated species GhTM-1 (G. hirsutum) and Gb3–79 (G. barbadense), and the fiber-less mutant GhMD17 (G. hirsutum background) were grown in soil in a greenhouse under natural sun light at The University of Texas at Austin. Flowers of cotton were tagged on the day of (post) anthesis (0 DPA) and collected at 1–3 pm to minimize diurnal effects. Immediately after flower collection, ovules were manually dissected from flowers under a dissection microscope. Flowers were collected in March 2016 for LCM RNA-seq analysis and in August 2018 for ovule culture assays. Unless noted otherwise, three biological replicates were used for each analysis.
Image analysis of cotton ovules and fiber cells
Cotton ovules were fixed with FAA solution (3.7% formaldehyde, 5% acetic acid, and 50% ethanol). They were then processed using an ethanol dilution series (50,70, 90, 95, and 100%), followed by an xylene:ethanol series (1:3, 1:1, 3:1, and 1:0), and finally multiple changes in a paraffin solution to embed siliques. Paraffin-embedded tissues were sectioned at 7 μm by microtome and mounted on glass side. Slides were deparaffinized using a brief xylene bath. Photos of ovules were taken using a light microscope ECLIPSE Ni (Nikon, Tokyo, Japan) and processed with ImageJ (https://imagej.nih.gov/ij/).
Laser-capture microdissection (LCM) and mRNA-seq library preparation
Ovules were collected from flowers at 0 and 2 DPA from GhTM-1, Gb3–79, and GhMD17, each with two biological replicates due to the laborious process, as performed in other studies using microarrays  and RNA-seq . Tissue fixation was performed as described in . Briefly, cotton ovaries were cut using a razor, and individual ovules were removed using forceps, fixed in a solution containing 75% ethanol and 25% acetic acid under vacuum for 15 min, and incubated at 4 °C for overnight. The samples transferred into fresh 75% ethanol and 25% acetic acid for 15 min on ice and microwaved at 250 W for 15 min using a PELCO BioWave (TED PELLA Inc., Redding, CA); the process was repeated three times each with a new solution. The samples were then dehydrated by ethanol series (70, 80, 90, 100, and 100%) with microwaving at 58 °C, 300 W, for 1.5 min at each step. Samples were then transferred to butanol in two steps (50%:50% ethanol:butanol followed by 100% butanol, each time microwaved at 58 °C, 300 W, for 1.5 min. Finally, samples were transferred into molten paraffin (50%:50% butanol:melted paraffin, then 100% molten paraffin, with microwaving at 58 °C, 250 W, for 10 min for each step). Molten paraffin was changed four times, with the sample vials immersed in a hot water bath, and microwaved at 58 °C, 250 W, for 30 min after each change. Finally, the samples were embedded in paraffin tissue blocks. Paraffin-embedded tissues were sectioned by a microtome (at 7 μm) and mounted on a pen-foil slide (ZEISS, Oberkochen, Germany). Fiber and epidermal cells were isolated and collected using Zeiss PALM Laser Microdissection with AdhesiveCap 500 (ZEISS, Oberkochen, Germany). The LCM fiber and epidermal tissues were subject to RNA extraction by RNA Aqueous Micro Total RNA Isolation Kit (Thermo Fisher Scientific, Waltham, MA), following the manufacturer’s protocol. The total RNA was treated with DNase treatment and qualified by a Bioanalyzer (Agilent, Santa Clara, CA) with RNA Integrity Number (RIN) > 6. An aliquot of total RNA (~ 100 ng) was used for RNA-seq library construction using the NEB Next Ultra RNA Library Prep Kit (New England Biolabs, Ipswich, MA). mRNA-seq libraries were sequenced by paired-end reads (150 bp) using an Illumina HiSeq 2500 (Illumina, San Diego, CA).
mRNA-seq data analysis of LCM samples
LCM mRNA-seq reads were mapped to G. hirsutum TM-1 genome v2.0 (GhTM-1 and GhMD17) and G. barbadense genome V1.0 (Gb3–79) whole genome  using Tophat 2.1.1 . Uniquely mapped reads were extracted and analyzed by Cufflinks 2.2.1  to determine gene expression levels using Fragments Per Kilobase Million (FPKM). Ortholog tables can be downloaded from Phytozome (https://genome.jgi.doe.gov/portal/pages/dynamicOrganismDownload.jsf?organism=Ghirsutum with the file name: /Ghirsutum/v2.1/orthology/Ghirsutum/inparanoid_Ghirsutum_527_v2.1.tar.gz/inParanoid_Ghirsutum_527_v2.1_Gbarbadense_526_v1.1.) (Supplemental Table 5a). Orthologs of G. hirsutum and G. barbadense (52,514 genes) expressed in these samples were used for further analysis. Values from two biological replicates were averaged to determine differentially expressed genes (DEGs) using the criteria of > 5 FPKM, ANOVA (P < 0.01), and fold changes > 1.5 (wild type epidermal vs. mutant epidermal; fiber vs. wild type epidermal) or > 2 (fiber vs. mutant epidermal).
Principal component analysis (PCA) was performed using log2-transformed gene expression values with singular value decomposition via the prcomp function in R (https://www.r-project.org).
Averaged gene expression values from two biological replicates were used for correlation coefficiency analysis. Pearsons’s correlation coefficients were determined using the cor function in R (https://www.r-project.org).
Gene Ontology (GO) analysis for upregulated and down-regulated genes was performed using the Bioconductor package topGO (https://bioconductor.org/packages/release/bioc/html/topGO.html) to identify overrepresented GO terms with statistical significance (P < 0.0001).
Expression values of all 15,220 homoeolog pairs were divided into seven groups as previously described : no expression (FPKM = 0 in both A and D), A = D (FPKM > 1 in both A and D, A < 2xD and D < 2xA), low expression (FPKM< 1 in both A and D), D = 0 (A FPKM> 1 and D = 0), A = 0 (D FPKM> 1 and A = 0), D > A (D is FPKM> 1 and D > 2xA), and A > D (A is FPKM> 1 and A > 2xD).
For phylogenetic analysis, protein sequences were obtained from Phyotozome using BIOMART, aligned using Clustal Omega  at EMBL-EBI (https://www.ebi.ac.uk/Tools/msa/clustalo/); the alignments were trimmed, and a neighbor-joining trees generated using JalView (http://www.jalview.org). Trees were visualized using Bioconductor package ggtree .
Ovules were collected from fertilized flowers at 0 DPA of GhTM-1 and GhMD17. Ovule culture was conducted using a modified protocol as previously described . Briefly, ovaries were isolated from flowers and surface sterilized in 100% ethanol for 15 min, and further dissected to expose the ovules. Individual ovules were collected and floated in BT medium  with or without 40 μg/L Rbin-1 (Axon Medchem, Groningen, Netherlands)  or 2 μg/L Roscovitine (#R7772-1MG, MilliporeSigma, Burlington, MA)  under dark conditions at 28 °C. For fiber initialization density counts, after 2 days, cultured ovules were briefly dried on a laboratory wipe before being gently placed on a microscope slide covered a UV curable adhesive (All-Ways Adhesives, Bensenville, IL). These slides were exposed to a hand-held UV torch for 30 s to set the adhesive before removing the ovules, and the resulting impressions were imaged via stereomicroscope SMZ1500 (Nikon, Tokyo, Japan). Fiber density was calculated by defining a region of interest in ImageJ, measuring the area and counting the fiber cell impressions.
For area measurements, cultured ovules were incubated for 21 days and then transferred to a glass plate where a stream of distilled water from a lab wash bottle was used to gently splay fibers. Excess water was then removed using laboratory wipes, and photos were taken. Total fiber and ovule areas were calculated and analyzed using ImageJ (https://imagej.nih.gov/ij/).
Availability of data and materials
The mRNA sequencing data have been deposited in the NCBI Gene Expression Omnibus (GEO) database, https://www.ncbi.nlm.nih.gov/geo (accession GSE152995).
Days post anthesis
TARGET OF RAPAMYCIN
Fragments per kilobase per million
Principal component analysis
ErbB-3 EPIDERMAL GROWTH FACTOR BINDING PROTEIN1
LETHAL WITH SEC THIRTEEN PROTEIN8
Chen ZJ, Scheffler BE, Dennis E, Triplett BA, Zhang T, Guo W, et al. Toward sequencing cotton (Gossypium) genomes. Plant Physiol. 2007;145(4):1303–10. https://doi.org/10.1104/pp.107.107672.
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.
Wendel JF, Cronn RC. Polyploidy and the evolutionary history of cotton. Adv Agron. 2003;78:139–86. https://doi.org/10.1016/S0065-2113(02)78004-8.
Chen ZJ, Sreedasyam A, Ando A, Song Q, De Santiago LM, Hulse-Kemp AM, et al. Genomic diversifications of five Gossypium allopolyploid species and their impact on cotton improvement. Nat Genet. 2020;52(5):525–33. https://doi.org/10.1038/s41588-020-0614-5.
Grover CE, Gallagher JP, Jareczek JJ, Page JT, Udall JA, Gore MA, et al. Re-evaluating the phylogeny of allopolyploid Gossypium L. Mol Phylogenet Evol. 2015;92:45–52. https://doi.org/10.1016/j.ympev.2015.05.023.
Splitstoser JC, Dillehay TD, Wouters J, Claro A. Early pre-Hispanic use of indigo blue in Peru. Sci Adv. 2016;2(9):e1501623. https://doi.org/10.1126/sciadv.1501623.
Kulkarni VN, Khadi BM, Maralappanavar MS, Deshapande LA, Narayanan SS. The Worldwide Gene Pools of Gossypium arboreum L. and G. herbaceum L., and Their Improvement. In: Paterson AH, editor. Genetics and Genomics of Cotton Plant Genetics and Genomics: Crops and Models, Volume 3, vol. 3. New York: Springer; 2009. p. 69–97.
Ruan YL. Boosting seed development as a new strategy to increase cotton fiber yield and quality. J Integr Plant Biol. 2013;55(7):572–5. https://doi.org/10.1111/jipb.12074.
Lee JJ, Hassan OS, Gao W, Wei NE, Kohel RJ, Chen XY, et al. Developmental and gene expression analyses of a cotton naked seed mutant. Planta. 2006;223(3):418–32. https://doi.org/10.1007/s00425-005-0098-7.
Chen ZJ, Guan X. Auxin boost for cotton. Nat Biotechnol. 2011;29(5):407–9. https://doi.org/10.1038/nbt.1858.
Haigler CH, Betancur L, Stiff MR, Tuttle JR. Cotton fiber: a powerful single-cell model for cell wall and cellulose research. Front Plant Sci. 2012;3:104.
Zhang M, Zheng X, Song S, Zeng Q, Hou L, Li D, et al. Spatiotemporal manipulation of auxin biosynthesis in cotton ovule epidermal cells enhances fiber yield and quality. Nat Biotechnol. 2011;29(5):453–8. https://doi.org/10.1038/nbt.1843.
Basra A, Malik CP. Development of the cotton fiber. Int Rev Cytol. 1984;89:65–113. https://doi.org/10.1016/S0074-7696(08)61300-5.
Taliercio E, Hendrix B, Stewart JM. DNA content and expression of genes related to cell cycling in developing Gossypium hirsutum (Malvaceae) fibers. Am J Bot. 2005;92(12):1942–7. https://doi.org/10.3732/ajb.92.12.1942.
De Veylder L, Beeckman T, Inze D. The ins and outs of the plant cell cycle. Nat Rev Mol Cell Biol. 2007;8(8):655–65. https://doi.org/10.1038/nrm2227.
Jones AR, Forero-Vargas M, Withers SP, Smith RS, Traas J, Dewitte W, et al. Cell-size dependent progression of the cell cycle creates homeostasis and flexibility of plant cell size. Nat Commun. 2017;8(1):15060. https://doi.org/10.1038/ncomms15060.
Barow M. Endopolyploidy in seed plants. Bioessays. 2006;28(3):271–81. https://doi.org/10.1002/bies.20371.
Van't Hof J. Increased nuclear DNA content in developing cotton fiber cells. Am J Bot. 1999;86(6):776–9. https://doi.org/10.2307/2656698.
Wu Y, Machado AC, White RG, Llewellyn DJ, Dennis ES. Expression profiling identifies genes expressed early during lint fibre initiation in cotton. Plant Cell Physiol. 2006;47(1):107–27. https://doi.org/10.1093/pcp/pci228.
Song Q, Guan X, Chen ZJ. Dynamic roles for small RNAs and DNA methylation during ovule and Fiber development in Allotetraploid cotton. PLoS Genet. 2015;11(12):e1005724. https://doi.org/10.1371/journal.pgen.1005724.
Pang M, Woodward AW, Agarwal V, Guan X, Ha M, Ramachandran V, et al. Genome-wide analysis reveals rapid and dynamic changes in miRNA and siRNA sequence and expression during ovule and fiber development in allotetraploid cotton (Gossypium hirsutum L). Genome Biol. 2009;10(11):R122. https://doi.org/10.1186/gb-2009-10-11-r122.
Wang M, Sun R, Li C, Wang Q, Zhang B. MicroRNA expression profiles during cotton (Gossypium hirsutum L) fiber early development. Sci Rep. 2017;7(1):44454. https://doi.org/10.1038/srep44454.
Donati G, Montanaro L, Derenzini M. Ribosome biogenesis and control of cell proliferation: p53 is not alone. Cancer Res. 2012;72(7):1602–7. https://doi.org/10.1158/0008-5472.CAN-11-3992.
Prakash V, Carson BB, Feenstra JM, Dass RA, Sekyrova P, Hoshino A, et al. Ribosome biogenesis during cell cycle arrest fuels EMT in development and disease. Nat Commun. 2019;10(1):2110. https://doi.org/10.1038/s41467-019-10100-8.
Lee JJ, Woodward AW, Chen ZJ. Gene expression changes and early events in cotton fibre development. Ann Bot (Lond). 2007;100(7):1391–401. https://doi.org/10.1093/aob/mcm232.
Xiao G, Zhao P, Zhang Y. A pivotal role of hormones in regulating cotton Fiber development. Front Plant Sci. 2019;10:87. https://doi.org/10.3389/fpls.2019.00087.
Wang M, Tu L, Yuan D, Zhu, Shen C, Li J, et al. Reference genome sequences of two cultivated allotetraploid cottons, Gossypium hirsutum and Gossypium barbadense. Nat Genet. 2019;51(2):224–9. https://doi.org/10.1038/s41588-018-0282-x.
Sun Y, Veerabomma S, Abdel-Mageed HA, Fokar M, Asami T, Yoshida S, et al. Brassinosteroid regulates fiber development on cultured cotton ovules. Plant Cell Physiol. 2005;46(8):1384–91. https://doi.org/10.1093/pcp/pci150.
Shi YH, Zhu SW, Mao XZ, Feng JX, Qin YM, Zhang L, et al. Transcriptome profiling, molecular biological, and physiological studies reveal a major role for ethylene in cotton fiber cell elongation. Plant Cell. 2006;18(3):651–64. https://doi.org/10.1105/tpc.105.040303.
Kim HJ, Hinchliffe DJ, Triplett BA, Chen ZJ, Stelly DM, Yeater KM, et al. Phytohormonal networks promote differentiation of fiber initials on pre-anthesis cotton ovules grown in vitro and in planta. PLoS One. 2015;10(4):e0125046. https://doi.org/10.1371/journal.pone.0125046.
Zeng J, Zhang M, Hou L, Bai W, Yan X, Hou N, et al. Cytokinin inhibits cotton fiber initiation by disrupting PIN3a-mediated asymmetric accumulation of auxin in the ovule epidermis. J Exp Bot. 2019;70(12):3139–51. https://doi.org/10.1093/jxb/erz162.
Davis LA, Addicott FT. Abscisic acid: correlations with abscission and with development in the cotton fruit. Plant Physiol. 1972;49(4):644–8. https://doi.org/10.1104/pp.49.4.644.
Marks MD. Molecular genetic analysis of trichome development in Arabidopsis. Ann Rev Plant Physiol Plant Mol Biol. 1997;48:137–63.
Guan X, Pang M, Nah G, Shi X, Ye W, Stelly DM, et al. miR828 and miR858 regulate homoeologous MYB2 gene functions in Arabidopsis trichome and cotton fibre development. Nat Commun. 2014;5(1):3050. https://doi.org/10.1038/ncomms4050.
Guan X, Lee JJ, Pang M, Shi X, Stelly DM, Chen ZJ. Activation of Arabidopsis seed hair development by cotton fiber-related genes. PLoS One. 2011;6(7):e21301. https://doi.org/10.1371/journal.pone.0021301.
Wang S, Wang JW, Yu N, Li CH, Luo B, Gou JY, et al. 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.
Walford SA, Wu Y, Llewellyn DJ, Dennis ES. GhMYB25-like: a key factor in early cotton fibre development. Plant J. 2011;65(5):785–97. https://doi.org/10.1111/j.1365-313X.2010.04464.x.
Zhang T, Hu Y, Jiang W, Fang L, Guan X, Chen J, et al. Sequencing of allotetraploid cotton (Gossypium hirsutum L. acc. TM-1) provides a resource for fiber improvement. Nat Biotechnol. 2015;33(5):531–7. https://doi.org/10.1038/nbt.3207.
Song Q, Zhang T, Stelly DM, Chen ZJ. Epigenomic and functional analyses reveal roles of epialleles in the loss of photoperiod sensitivity during domestication of allotetraploid cottons. Genome Biol. 2017;18(1):99. https://doi.org/10.1186/s13059-017-1229-8.
Turley RB, Kloth RH. Identification of a third fuzzless seed locus in upland cotton (Gossypium hirsutum L.). J Hered. 2002;93(5):359–64. https://doi.org/10.1093/jhered/93.5.359.
Wan Q, Guan X, Yang N, Wu H, Pan M, Liu B, et al. Small interfering RNAs from bidirectional transcripts of GhMML3_A12 regulate cotton fiber development. New Phytol. 2016;210(4):1298–310. https://doi.org/10.1111/nph.13860.
Zhu QH, Yuan Y, Stiller W, Jia Y, Wang P, Pan Z, et al. Genetic dissection of the fuzzless seed trait in Gossypium barbadense. J Exp Bot. 2018;69(5):997–1009. https://doi.org/10.1093/jxb/erx459.
Clark G, Torres J, Finlayson S, Guan X, Handley C, Lee J, et al. Apyrase (nucleoside triphosphate-diphosphohydrolase) and extracellular nucleotides regulate cotton fiber elongation in cultured ovules. Plant Physiol. 2010;152(2):1073–83. https://doi.org/10.1104/pp.109.147637.
Saez-Vasquez J, Delseny M. Ribosome biogenesis in plants: from functional 45S ribosomal DNA organization to ribosome assembly factors. Plant Cell. 2019;31(9):1945–67. https://doi.org/10.1105/tpc.18.00874.
Xue S, Barna M. Specialized ribosomes: a new frontier in gene regulation and organismal biology. Nat Rev Mol Cell Biol. 2012;13(6):355–69. https://doi.org/10.1038/nrm3359.
Ferreyra ML, Biarc J, Burlingame AL, Casati P. Arabidopsis L10 ribosomal proteins in UV-B responses. Plant Signal Behav. 2010;5(10):1222–5. https://doi.org/10.4161/psb.5.10.12758.
Williams ME, Sussex IM. Developmental regulation of ribosomal protein L16 genes in Arabidopsis thaliana. Plant J. 1995;8(1):65–76. https://doi.org/10.1046/j.1365-313X.1995.08010065.x.
Weijers D, Franke-van Dijk M, Vencken RJ, Quint A, Hooykaas P, Offringa R. An Arabidopsis minute-like phenotype caused by a semi-dominant mutation in a RIBOSOMAL PROTEIN S5 gene. Development. 2001;128(21):4289–99.
Wullschleger S, Loewith R, Hall MN. TOR signaling in growth and metabolism. Cell. 2006;124(3):471–84. https://doi.org/10.1016/j.cell.2006.01.016.
Blenis J. TOR, the gateway to cellular metabolism, cell growth, and disease. Cell. 2017;171(1):10–3. https://doi.org/10.1016/j.cell.2017.08.019.
Shi L, Wu Y, Sheen J. TOR signaling in plants: conservation and innovation. Development. 2018;145(13):dev160887.
Menand B, Desnos T, Nussaume L, Berger F, Bouchez D, Meyer C, et al. Expression and disruption of the Arabidopsis TOR (target of rapamycin) gene. Proc Natl Acad Sci U S A. 2002;99(9):6422–7. https://doi.org/10.1073/pnas.092141899.
Horvath BM, Magyar Z, Zhang Y, Hamburger AW, Bako L, Visser RG, Bachem CW, Bogre L, EBP1 regulates organ size through cell growth and proliferation in plants. EMBO J. 2006;25(20):4909–20. https://doi.org/10.1038/sj.emboj.7601362.
Kawashima SA, Chen Z, Aoi Y, Patgiri A, Kobayashi Y, Nurse P, et al. Potent, reversible, and specific chemical inhibitors of eukaryotic ribosome biogenesis. Cell. 2016;167(2):512–24 e514. https://doi.org/10.1016/j.cell.2016.08.070.
Yarotskyy V, Elmslie KS. Roscovitine, a cyclin-dependent kinase inhibitor, affects several gating mechanisms to inhibit cardiac L-type (ca(V)1.2) calcium channels. Br J Pharmacol. 2007;152(3):386–95. https://doi.org/10.1038/sj.bjp.0707414.
Beasley CA, Ting IP. The effects of plant growth substances on in vitro fiber development from unfertilized cotton ovules. Amer J Bot. 1974;61(2):188–94. https://doi.org/10.1002/j.1537-2197.1974.tb06045.x.
Gilbert MK, Bland JM, Shockey JM, Cao H, Hinchliffe DJ, Fang DD, et al. A transcript profiling approach reveals an abscisic acid-specific glycosyltransferase (UGT73C14) induced in developing fiber of Ligon lintless-2 mutant of cotton (Gossypium hirsutum L.). PLoS One. 2013;8(9):e75268.
Saez A, Apostolova N, Gonzalez-Guzman M, Gonzalez-Garcia MP, Nicolas C, Lorenzo O, et al. Gain-of-function and loss-of-function phenotypes of the protein phosphatase 2C HAB1 reveal its role as a negative regulator of abscisic acid signalling. Plant J. 2004;37(3):354–69. https://doi.org/10.1046/j.1365-313X.2003.01966.x.
Dupeux F, Santiago J, Betz K, Twycross J, Park SY, Rodriguez L, et al. A thermodynamic switch modulates abscisic acid receptor sensitivity. EMBO J. 2011;30(20):4171–84. https://doi.org/10.1038/emboj.2011.294.
Zhang Y, He P, Yang Z, Huang G, Wang L, Pang C, et al. A genome-scale analysis of the PIN gene family reveals its functions in cotton Fiber development. Front Plant Sci. 2017;8:461.
Xiao YH, Li DM, Yin MH, Li XB, Zhang M, Wang YJ, et al. Gibberellin 20-oxidase promotes initiation and elongation of cotton fibers by regulating gibberellin synthesis. J Plant Physiol. 2010;167(10):829–37. https://doi.org/10.1016/j.jplph.2010.01.003.
Guo H, Ecker JR. The ethylene signaling pathway: new insights. Curr Opin Plant Biol. 2004;7(1):40–9. https://doi.org/10.1016/j.pbi.2003.11.011.
Song Q, Ando A, Xu D, Fang L, Zhang T, Huq E, et al. Diurnal down-regulation of ethylene biosynthesis mediates biomass heterosis. Proc Natl Acad Sci U S A. 2018;115(21):5606–11. https://doi.org/10.1073/pnas.1722068115.
Luo M, Xiao Y, Li X, Lu X, Deng W, Li D, et al. GhDET2, a steroid 5alpha-reductase, plays an important role in cotton fiber cell initiation and elongation. Plant J. 2007;51(3):419–30. https://doi.org/10.1111/j.1365-313X.2007.03144.x.
Zhou Y, Zhang ZT, Li M, Wei XZ, Li XJ, Li BY, et al. Cotton (Gossypium hirsutum) 14-3-3 proteins participate in regulation of fibre initiation and elongation by modulating brassinosteroid signalling. Plant Biotechnol J. 2015;13(2):269–80. https://doi.org/10.1111/pbi.12275.
Mora-Garcia S, Vert G, Yin Y, Cano-Delgado A, Cheong H, Chory J. Nuclear protein phosphatases with Kelch-repeat domains modulate the response to brassinosteroids in Arabidopsis. Genes Dev. 2004;18(4):448–60. https://doi.org/10.1101/gad.1174204.
Yin Y, Wang ZY, Mora-Garcia S, Li J, Yoshida S, Asami T, et al. BES1 accumulates in the nucleus in response to brassinosteroids to regulate gene expression and promote stem elongation. Cell. 2002;109(2):181–91. https://doi.org/10.1016/S0092-8674(02)00721-3.
Machado A, Wu Y, Yang Y, Llewellyn DJ, Dennis ES. The MYB transcription factor GhMYB25 regulates early fibre and trichome development. Plant J. 2009;59(1):52–62. https://doi.org/10.1111/j.1365-313X.2009.03847.x.
Loewith R, Hall MN. Target of rapamycin (TOR) in nutrient signaling and growth control. Genetics. 2011;189(4):1177–201. https://doi.org/10.1534/genetics.111.133363.
Saha S, Hof JVT. Cotton Fiber cells are arrested at G1 stage. J New Seeds. 2005;7(1):1–8. https://doi.org/10.1300/J153v07n01_01.
Xiao G, He P, Zhao P, Liu H, Zhang L, Pang C, et al. Genome-wide identification of the GhARF gene family reveals that GhARF2 and GhARF18 are involved in cotton fibre cell initiation. J Exp Bot. 2018;69(18):4323–37. https://doi.org/10.1093/jxb/ery219.
Kravchenko A, Citerne S, Jehanno I, Bersimbaev RI, Veit B, Meyer C, et al. Mutations in the Arabidopsis Lst8 and Raptor genes encoding partners of the TOR complex, or inhibition of TOR activity decrease abscisic acid (ABA) synthesis. Biochem Biophys Res Commun. 2015;467(4):992–7. https://doi.org/10.1016/j.bbrc.2015.10.028.
Schepetilnikov M, Dimitrova M, Mancera-Martinez E, Geldreich A, Keller M, Ryabova LA. TOR and S6K1 promote translation reinitiation of uORF-containing mRNAs via phosphorylation of eIF3h. EMBO J. 2013;32(8):1087–102. https://doi.org/10.1038/emboj.2013.61.
Schepetilnikov M, Ryabova LA. Auxin signaling in regulation of plant translation Reinitiation. Front Plant Sci. 2017;8:1014. https://doi.org/10.3389/fpls.2017.01014.
Oshima Y, Shikata M, Koyama T, Ohtsubo N, Mitsuda N, Ohme-Takagi M. MIXTA-like transcription factors and WAX INDUCER1/SHINE1 coordinately regulate cuticle development in Arabidopsis and Torenia fournieri. Plant Cell. 2013;25(5):1609–24. https://doi.org/10.1105/tpc.113.110783.
Gilding EK, Marks MD. Analysis of purified glabra3-shapeshifter trichomes reveals a role for NOECK in regulating early trichome morphogenic events. Plant J. 2010;64(2):304–17. https://doi.org/10.1111/j.1365-313X.2010.04329.x.
Belmonte MF, Kirkbride RC, Stone SL, Pelletier JM, Bui AQ, Yeung EC, et al. Comprehensive developmental profiles of gene activity in regions and subregions of the Arabidopsis seed. Proc Natl Acad Sci U S A. 2013;110(5):E435–44. https://doi.org/10.1073/pnas.1222061110.
Henry KF, Goldberg RB. Using giant scarlet runner bean embryos to uncover regulatory networks controlling suspensor gene activity. Front Plant Sci. 2015;6. https://doi.org/10.3389/fpls.2015.00044.
Takahashi H, Kamakura H, Sato Y, Shiono K, Abiko T, Tsutsumi N, et al. A method for obtaining high quality RNA from paraffin sections of plant tissues by laser microdissection. J Plant Res. 2010;123(6):807–13. https://doi.org/10.1007/s10265-010-0319-4.
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(4).
Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, et al. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and cufflinks. Nat Protoc. 2012;7(3):562–78. https://doi.org/10.1038/nprot.2012.016.
Sievers F, Wilm A, Dineen D, Gibson TJ, Karplus K, Li W, et al. Fast, scalable generation of high-quality protein multiple sequence alignments using Clustal omega. Mol Syst Biol. 2011;7(1):539. https://doi.org/10.1038/msb.2011.75.
Yu G, Smith DK, Zhu H, Guan Y, Lam TT. GGTREE: an R package for visualization and annotation of phylogenetic trees with their covariates and otherassociated data. Methods Ecol Evol. 2017;8(1):28–36. https://doi.org/10.1111/2041-210X.12628.
We thank Mr. Shane Merrell of Plant Growth Facilities and Ms. Julie Hayes and Dr. Dwight Romanovicz of Center for Biomedical Research Support at The University of Texas at Austin for assistance in plant growth, microscopic and imaging analysis and laser-capture microdissection, and Texas Advanced Computing Center for computational support and assistance.
This work is supported by grants from the National Science Foundation (IOS1739092) and the Cotton Incorporated (20–799).
Ethics approval and consent to participate
Consent for publication
The authors declare no conflict of interest in this work.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Images of Laser Capture Microdissection (LCM) and RNA quality in LCM samples. a Isolation of fiber and epidermal cell tissues by Laser Capture Microdissection (LCM) method. Top represents epidermal cell isolation from ovules at 0 DPA in GhTM-1 and Gb3–79, and GhMD17. Bottom shows separation of fiber and epidermal cell layers from ovules at 2 DPA in GhTM-1 and Gb3–79. Bar represents 150 μm. b Line plots showing RNA size (x-axis) and quantity (y-axis) measured by Bioanalyzer with RIN number in GhTM-1 epidermal cells (0 DPA) (left) and GhTM-1 fiber cells (2 DPA) (right). Two ribosomal RNA peaks are clearly visible. Supplemental Figure 2. Gene expression correlation between GhTM-1, Gb3–79, and GhMD17 and down-regulated genes in early fiber cells. a Pearson’s correlation coefficients of LCM RNA-seq samples averaged between 2 replicates. Hierarchical clustering was based on geometric distance. b Venn diagram analysis of down-regulated genes in early fiber cells at 0 DPA; GhMD17 > GhTM-1 (blue) and GhMD17 > Gb3–79 (red). Shared down-regulated genes in fibers at 0 DPA (white) between comparisons are shown in dashed box. c Venn diagram analysis of down-regulated genes in early fiber cells at 2 DPA; GhTM-1(epi > fiber) (blue), Gh3–79 (epi > fiber) (red), GhMD17(epi) > GhTM-1(fiber) (yellow), and GhMD17(epi) > Gb3–79(fiber) (green). Shared down-regulated genes in early fiber cells at 2 DPA (white) between comparisons are shown in dashed box. Supplemental Figure 3. Gene expression divergence between GhTM-1 and Gb3–79 in epidermal cells at 0 and 2 DPA. a Venn Diagram analysis of upregulated genes in epidermal cells at 0 DPA; GhTM-1 0epi > 2epi (light blue) and Gb3–79 0epi > 2epi (pink). b Venn Diagram analysis of upregulated genes in epidermal cells at 2 DPA; GhTM-1 2epi > 0epi (blue) and Gb3–79 2epi > 0epi (red). Supplemental Figure 4. Biased expression of A and D homoeologs in GhTM-1 and Gb3–79. a Number of A and D homoeologs (pairs), showing no expression (orange), A = D (red), low expression (FPKM< 1) (yellow), D = 0 (A FPKM> 1) (light green), A = 0 (D FPKM> 1) (green), D > A (D FPKM> 1 and D > 2xA) (blue), and A > D (A FPKM> 1 and A > 2xD) (light blue). b Number of A and D homoeologs (pairs), showing D = 0 (A FPKM> 1) (light green), A = 0 (D FPKM> 1) (green), D > A (D FPKM> 1 and D > 2xA) (blue), and A > D (A FPKM> 1 and A > 2xD) (light blue). c Number of A and D homoeologs with biased expression, either A-homoeologs highly expressed (red) or D-homoeologs highly expressed (blue). d Percentage of A and D homoeologs (pairs) with biased expression involved in cell cycle regulation (251 homeolog pairs) and ribosomal protein biosynthesis (81 homeolog pairs) between epidermal (0 and 2 DPA) and fiber cells. Colors represent A homoeologs (cell cycle, red; ribosome, pink) and D homoeologs (cell cycle, blue; ribosome, light blue). Supplemental Figure 5. Clustering analysis of expressed genes encoding cell cycle regulators and ribosomal protein subunits. a Expression clustering of cell cycle related genes in 0 DPA epidermal (left), 2 DPA epidermal (middle), 2 DPA fiber (right) of GhTM-1, Gb3–79, and GhMD17. Upregulation (yellow box) and down-regulation (green box) of genes in early fiber cells are indicated. b Expression clustering of putative ribosome protein subunit genes in 0 DPA epidermal (left), 2 DPA epidermal (middle), 2 DPA fiber (right) of GhTM-1, Gb3–79, and GhMD17. Upregulation (yellow box) and down-regulation (green box) in early fiber cells are indicated. Supplemental Figure 6. Analysis of fiber cell numbers in cultured ovules. a Number of fiber cells per 1 mm2 on ovule surface of GhTM-1 incubated in culture for 2 days with control (grey), Rbin-1 (blue), and Roscovitine (red) treatments. Two asterisks indicate a statistical significance level of P < 0.01 (Student’s t test). b Photos of epidermis showing fiber cell initials in GhTM-1 incubated with mock (control, top), Rbin-1 (middle), and Roscovitine (bottom). Bars = 60 μm. Supplemental Figure 7. Phylogenetic analysis of upregulated and down-regulated ARFs in cotton fibers and their Arabidopsis homologs. Phylogenetic tree of all ARF genes in cotton and Arabidopsis. Gene names in the figure were identified in the list of upregulated and down-regulated genes in cotton fibers and their homologs in Arabidopsis. Color of line in phylogenetic tress represents upregulated (red), down-regulated (blue) genes, both up and down-regulated genes (green, ARF19) in fibers. Supplemental Figure 8. Phylogenetic analysis of upregulated and down-regulated MYB transcription factor genes in cotton fibers and their Arabidopsis homologs. Phylogenetic trees of all MYB transcription factor genes in cotton and Arabidopsis. Gene names in the figure were identified in the list of upregulated and down-regulated genes in cotton early fiber cells and their homologs in Arabidopsis. Genes described and named in previous Gossypium studies listed in parenthesis. Color of line in phylogenetic tress represents up-regulated (red) and down-regulated (blue) genes in fibers, green are other genes known to be upregulated in fibers which were not detected in this experiment.
RNA-seq statistics of LCM samples. Supplemental Table 2. Upregulated and down-regulated genes in fiber cell initials at 0 DPA. Supplemental Table 3. Upregulated and down-regulated genes in early fiber cells at 2 DPA. Supplemental Table 3a. Up-regulated and down-regulated genes in epidermal cells in 0 and 2 DPA. Supplemental Table 4. Upregulated genes in early fiber cells (0 and 2 DPA) specific to species. Supplemental Table 4a. Upregulated and down-regulated genes in fiber cells between GhTM-1 and Gb3–79 at 2 DPA. Supplemental Table 4b. Significant GO terms in upregulated and down-regulated genes between GhTM-1 and Gb3–79 in fiber cells at 2 DPA. Supplemental Table 5. List of A and D homoeologs with biased expression values. Supplemental Table 5a. Orthologs of G. hirstum and G. barbadense. Supplemental Table 6. List of putative ribosomal protein subunit and cell cycle genes with expression values in LCM samples. Supplemental Table 7. Expression values of putative TOR, RAPTOR, LST8, and EPB1 genes in LCM samples. Supplemental Table 8. Expression values of phytohormone and MYB transcription factor related genes in LCM samples.
About this article
Cite this article
Ando, A., Kirkbride, R.C., Jones, D.C. et al. LCM and RNA-seq analyses revealed roles of cell cycle and translational regulation and homoeolog expression bias in cotton fiber cell initiation. BMC Genomics 22, 309 (2021). https://doi.org/10.1186/s12864-021-07579-1