Novel insights into water-deficit-responsive mRNAs and lncRNAs during fiber development in Gossypium hirsutum

Abstract

Background: The fiber yield and quality of cotton are greatly and periodically affected by water deficit. However, the molecular mechanism of the water deficit response in cotton fiber cells has not been fully elucidated.

Results: In this study, water deficit caused a significant reduction in fiber length, strength, and elongation rate but a dramatic increase in micronaire value. To explore genome-wide transcriptional changes, fibers from cotton plants subjected to water deficit (WD) and normal irrigation (NI) during fiber development were analyzed by transcriptome sequencing. Analysis showed that 3427 mRNAs and 1021 long noncoding RNAs (lncRNAs) from fibers were differentially expressed between WD and NI plants. The maximum number of differentially expressed genes (DEGs) and lncRNAs (DERs) was identified in fibers at the secondary cell wall biosynthesis stage, suggesting that this is a critical period in response to water deficit. Twelve genes in cotton fiber were differentially and persistently expressed at ≥ five time points, suggesting that these genes are involved in both fiber development and the water-deficit response and could potentially be used in breeding to improve cotton resistance to drought stress. A total of 540 DEGs were predicted to be potentially regulated by DERs by analysis of coexpression and genomic colocation, accounting for approximately 15.76% of all DEGs. Four DERs, potentially acting as target mimics for microRNAs (miRNAs), indirectly regulated their corresponding DEGs in response to water deficit.

Conclusions: This work provides a comprehensive transcriptome analysis of fiber cells and a set of protein-coding genes and lncRNAs implicated in the cotton response to water deficit, significantly affecting fiber quality during the fiber development stage.

Keywords: Fiber; Gossypium hirsutum; Long noncoding RNA (lncRNA); Water deficit; mRNA.

Fig. 1

Water deficit caused a significant reduction in fiber quality. A SRWC in the NI-field and WD-field at 0, 15, 30 and 35 DPA. Seed cotton weight (B), seed index (C), lint weight (D), lint index (E), lint (F), length uniformity (G), fiber length (H), fiber strength (I), fiber elongation rate (J), and micronaire (K), for NI-treatments were compared with WD-treatments. Data represent the mean ± SE of seven biological replicates. Fibers for each replication were sampled from 20 naturally-open bolls on the middle section of cotton plants. Confidence levels were tested by unpaired t-test (*, P < 0.05; ns, not significant)

Fig. 2

Identification and characterization of mRNAs and lncRNAs in G. hirsutum fibers. A Representative images of individual seeds with attached fibers from 0 DPA to 35 DPA. Cotton fibers undergo three major sequential and overlapping developmental stages before maturity, including initiation, elongation and secondary cell wall biosynthesis. Transition-1 and -2 are two fiber developmental transition stages, which are from initiation to elongation and from elongation to secondary cell wall biosynthesis, respectively. The scale bars in all panels are 0.5 cm. B The pipeline of mRNAs and lncRNAs identification. C Length density distributions of lincRNAs, lncNATs and mRNAs. D Exon number per transcript of lincRNAs, lncNATs and mRNAs. E FPKM distributions of lincRNAs, lncNATs and mRNAs

Fig. 3

Confirmation of the expression patterns of mRNAs (A) and lncRNAs (B) using qPCR. Ten mRNAs and ten lncRNAs were randomly selected for expression analysis during the fiber developmental stages of G. hirsutum ND13 treated with NI and WD. The correlation of relative expression for mRNAs and lncRNAs measured by RNA-seq and qPCR was estimated with r-values. UBQ7 was used as the reference gene. Gene (mRNA) IDs are shown in the genome of G. hirsutum TM-1 (NAU-NBI v1.1), including Gh_A05G0770 (17.3 kDa class I heat shock protein), Gh_A09G1977 (1-aminocyclopropane-1-carboxylate oxidase homolog 1), Gh_A11G2903 (ABC transporter G family member 2), Gh_D01G0047 (Protein RADIALIS-like 6), Gh_D03G1452 (Tubulin beta-7 chain), Gh_D04G0942 (No annotation), Gh_D05G1621 (No annotation), Gh_D08G1970 (Probable aquaporin PIP1-2), Gh_D08G2730 (Bidirectional sugar transporter SWEET15), and Gh_Sca115726G01 (Aspartic proteinase nepenthesin 1)

Fig. 4

Identification and characterization of DEGs between NI-treated and WD-treated cotton fibers. Twelve genes were differentially expressed at ≥5 time points, including ADH (alcohol dehydrogenase), MIOX (myo-inositol oxygenase), TK (thymidine kinase), PS (phosphate starvation-induced gene), PIP (plasmamembrane intrinsic protein), PAP (purple acid phosphatase), SPX (SYG1-Pho81-XPR1 domain-containing protein), NAM (no apical meristem), NCED (9-cis-epoxycarotenoid dioxygenase), UMAMIT (usually multiple acids move in and out transporter) and PEPC (phosphoethanolamine/phosphocholine phosphatase)

Fig. 5

The significantly enriched KEGG pathways of DEGs between NI-treated and WD-treated cotton fibers. The overall trends of upregulation and downregulation for DEGs are indicated by red and green arrows, respectively

Fig. 6

Identification and characterization of DERs between NI-treated and WD-treated cotton fibers. Green bars for downregulated DERs. Purple bars for upregulated DERs. Red bars and dots for DERs that were only shown at one time point

Fig. 7

DEGs-R predicted by gene coexpression and genomic colocation analysis for DEGs and DERs. Blue bars for DEGs. Orange bars for DEGs-R. Green dots for the percentages of DEGs-R/DEGs

Fig. 8

Prediction of the interactions between lncRNAs, miRNAs and mRNAs by forming RNA-RNA duplexes. miRNA-directed target mRNA degradation was potentially regulated by forming a lncRNA-miRNA duplex. Gene (mRNA) IDs are shown in the genome of G. hirsutum TM-1 (NAU-NBI v1.1), including Gh_A10G1972 (DEAD-box ATP-dependent RNA helicase 42), Gh_D06G2174 (protein of unknown function), Gh_A07G2348 and Gh_D07G0162 (LRR receptor-like serine/threonine-protein kinase), and Gh_A07G2019 (UDP-glycosyltransferase 88F3). The expression of lncRNAs and mRNAs is shown with log₂FoldChange (WD/NI)