The Cytoplasmic Male Sterility Biology Essay

Published: November 2, 2015 Words: 4819

The cytoplasmic male sterility (CMS), a system in which a plant is unable to produce fertile pollen, is an important tool for hybrid seed production without the need of hand emasculation to utilize heterosis in many crops. It is also an ideal biological system for studying male gametogenesis and nucleo-cytoplasmic interactions. CMS, as a maternally inherited trait, is associated with novel open reading frames (orfs) found in plant mitochondrial genomes. Therefore, dysfunctions of mitochondria during male meiosis and pollen development have been the focus to the understanding of the molecular mechanism of CMS. However, it remains unclear how these dysfunctions are associated with microRNAs (miRNAs).

Results

In this study, a comparative deep RNA sequencing was performed to identify differentially expressed miRNAs from cotton flower buds at the stage of male meiosis between CMS-D2 and its isogenic maintainer line with a fertile cytoplasm. Of more than 9-10 million short sequencing reads generated for each of the two lines, 43 conserved miRNAs were recovered, many of which have not been previously observed in cotton. In addition, 86 non-conserved novel miRNA candidates were predicted. Through a statistical analysis, 27 differentially expressed miRNAs were identified between the CMS-D2 and its maintainer line, indicating possible associations between CMS and miRNA regulations. Real-time PCR analyses further confirmed the differential expression profiles of miR156 and its target SQUAMOSA Promoter-Binding Protein-Like (SPL-like) genes between the two genotypes.

Conclusions

Affected by the male sterile cytoplasm, the expression of many miRNAs and their target genes was changed in CMS-D2. These miRNAs and their target genes may be involved in mitochondria-nuclear interactions. Thus, our study provides the first glimpse into miRNAs related to male meiosis and CMS for further studies of mitochondria-nuclear interactions in Upland cotton.

Introduction

Cotton is one of the most important economic crops in the world, and Upland cotton (Gossypium hirsutum) predominates the world's cotton fiber production. Several strategies have been developed for a better utilization of heterosis in cotton. The cytoplasmic male sterility (CMS) and restoration system is considered an important approach for hybrid cotton seed production, and also an ideal tool for studying male gametogenesis and nucleo-cytoplasmic interactions.

In Upland cotton, CMS-D2 with G. harknessii (D2-2 genome) cytoplasm and CMS-D8 with G. trilobum (D8 genome) cytoplasm are the two main CMS systems [1,2,3]. To identify the CMS-associated genes, Feng et al. [4] first reported a comparative restriction fragment length ploymorphism (RFLP) analysis of mitochondrial DNA among CMS-D2, its fertile restorer and maintainer. Recently, Wang et al. [5] and Wu et al. [6] demonstrated that mitochondrial DNA gene atpA in the CMS mitochondrial genomes may be one of the candidates for CMS-associated genes in CMS cotton. However, little is known regarding interactions and regulations between atpA and nuclear genes. To date, most studies in the area of CMS are focused on identifying the structural and functional changes in the mitochondria genome that cause male sterility in different plants [7]. While male meiosis and pollen development is mostly controlled by nuclear genes, little is known about how mitochondrial dysfunctions affect the reproductive development in a CMS system. As such, for a better understanding of the nucleo-cytoplasmic interaction in CMS, it is necessary to gain further insight into the downstream regulations under the CMS cytoplasm.

In plants, known miRNAs are produced as 21-22 nt duplexes from their precursors by Dicer-like 1 (DCL1), and they function in post-transcriptional gene silencing by guiding mRNA degradation or translational repression through perfect or near-perfect complementarity with target mRNAs [8,9,10]. Recently, 24-nt miRNAs produced by DCL3/RDR2/Pol IV pathway that participated in DNA methylation were also identified [11,12]. Now, more than 1000 miRNAs have been annotated in Arabidopsis, rice and other plant species [13]. Several experimental and genetic analyses have indicated that miRNAs play a key regulatory role in plant growth and development [14]. Most recently, evidence indicates that miRNAs are also present during pollen development [15,16], and play a vital role in gene regulation at various stages of pollen development [17]. However, it is unknown whether the expression profiles of miRNAs are affected by a CMS cytoplasm.

Considering the importance of miRNAs in regulating plant growth and development, a large number of putative miRNAs was identified in cotton by computational approaches and high throughput sequencing technology, and most of the reports were focused on ovule and fiber development [18,19,20,21]. The objectives of the current study were to identify differentially expressed miRNAs from cotton flower buds at the stage of male meiosis between CMS-D2 and its isogenic maintainer line, and to predict their possible target genes. The results represent one of the first reports in identification of miRNAs with regarding to CMS in plants.

Results

Overall quality of small RNA libraries

Using the Illumina sequencing technology, a total of 11,449,151 and 9,823,976 sequencing reads were generated, consisting of 5,202,242 and 4,150,875 sequences for the maintainer (ZB) and its CMS (ZBA) genotypes, respectively. Prior to mapping these sequences to reference sequences, a series of filtering were performed to ensure high quality sequences, including: (1) adaptor sequences filter (0.92% and 0.76% of raw data counts in ZB and ZBA, corresponding to 0.37% and 0.30% sequencing reads respectively), (2) low quality sequence filter (0.05% and 0.06% of sequencing reads in ZB and ZBA, corresponding to 0.04% and 0.05% of sequences, respectively), (3) sequence length filter (5.24% and 4.57% of sequencing reads in ZB and ZBA, corresponding to 5.72% and 5.11% of sequences, respectively), (4) simple sequence filter (0.42% and 0.36% of sequencing reads in ZB and ZBA, corresponding to 0.80% and 0.74% of sequences, respectively), and (5) sequence of copy number filter (41.67% and 39.01% of sequencing reads in ZB and ZBA, corresponding to 85.41% and 85.52% of sequences, respectively) (Supplement Table 1). In addition, certain known types of RNA sequences (e.g., mRNA, rRNA, tRNA, snRNA, snoRNA and repetitive sequence elements) were also excluded from analysis (Supplement Table 2). However, these sequences accounted for only a small fraction (1.33% and 1.56% for ZB and ZBA, respectively) of the sequences, but 17.63% and 23.24% in sequencing reads respectively. The detection of a low proportion of these RNAs indicated that the sequencing samples were not contaminated by degraded RNA, and were therefore of high integrity.

Of the mappable sequences in the two genotypes, the majority were 21-24 nt in length, which is typical of Dicer-processed small RNA (Fig. 1). In both genotypes, 21 nt small RNAs were the most abundant, followed by the 22 and 24 nt long small RNAs. It appeared that the total number of mappable sequences in the maintainer line was more than that in the CMS line.

Known cotton miRNAs

For the two lines, the mappable sequences were divided into three categories based on matching identity with the known miRNAs reported in cotton and other plant species, miRNA*s (partner strand of mature miRNA) and pre-miRNAs. All these sequences could find their homologous sequences in cotton genome survey sequence (GSS) or expressed sequence tag (EST) sequences in the NCBI database. Among the 33 known cotton miRNAs reported in miRBase database (http://www.mirbase.org/index.shtml, released 15.0), 11 miRNAs corresponding to 9 pre-miRNAs were investigated in this study. One cotton mature miRNA, ghr-miR167, which was cloned by Pang [21], was found to be expressed in the two genotypes, but we failed to find its pre-miRNA sequence because of lacking of sequence information in the GSS and EST database. The mature miRNAs typically have relatively higher steady-state levels than their corresponding miRNA*. However, some miRNA* known to have regulatory roles [22], reached substantial levels. In this study, two miRNA*s (PN -ghr-miR156a-3p and PN-ghr-miR164-3p, here PN means Predicted Novel) were not reported in cotton before (Supplement Table 3).

Conserved plant miRNAs in cotton

Supplement Table 4 lists a total of 32 miRNAs which were originally identified in other plants as pre-miRNAs but were unable to be mapped to cotton GSS and EST in the current study. Five miRNAs* (PC-ghr-miR167-3p, PC-ghr-miR156-3p, PC-ghr-miR319-3p, PC-ghr-miR482-3p and PC-ghr-miR168-3p; here PC means Plant Conserved) were first reported in Upland cotton.

Novel miRNA predicted

A total of 86 predicted novel (PN type) miRNAs corresponding to 86 pre-miRNAs were identified in the two genotypes (Supplement Table 5). Most (80.2%) of them were 24 nt in length. 23 predicted miRNAs were identified in the EST database and 63 were present in the GSS database. Diversity in locations of mature miRNA sequences was observed for the newly identified cotton PN-miRNAs in that 47 mature miRNA sequences are located at the 5' end of the PN-miRNA precursor sequences and others are at the 3' end. As with other plants, cotton miRNA precursors showed diversity in size and structure. For example, the length of miRNA precursors varied from 67 to 256 nt with an average of 139 nt. The sequences of the miRNA precursors have AU content ranging from 35.9 to 80% with an average of 62.9%, and 78 of 86 miRNA precursors have AU content higher than 50%. The results are in agreement with the notion that miRNA precursors and mature miRNAs contain more AU nucleotides than GC, enabling miRNA precursors to be more easily processed into mature miRNA by the RISC complex [23]. MFEI (minimal folding free energy index) is a useful criterion for distinguishing miRNAs from other types of coding and non-coding RNAs [23]. The newly identified cotton pre-miRNAs had MFEI values ranging from 0.7 to 1.3, with an average of 0.86. This satisfied the criterion that a sequence with MFEI > 0.85 is likely a miRNA.

Differentially expressed miRNAs and their predicted targeted genes

For sequence-based differential expression profiles, the total number of the mappable reads of all the isoforms of the miRNAs for each genotype was used as the scaling factor for identification of differential expressed (DE) miRNAs. A total of 27 DE miRNAs (5 known cotton miRNAs, 17 PC-type miRNAs and 5 PN-type miRNAs) were identified between the two genotypes using IDEG6 with a Chi-square test [24]. Eleven miRNAs were down-regulated and 16 miRNAs were up-regulated in the sterile CMS-D2 line, including three CMS specifically expressed miRNAs (PC-ghr-miR319, PN-ghr-60 and PN-ghr-72) and two fertile cytoplasm specifically expressed miRNAs (PN-ghr-64 and PN-ghr-66) (Table 1).

Previous studies have shown that miRNAs regulate gene expression by binding to targeted mRNA sequences in a perfect or ''near-perfect'' complementary site in plants [25]. This makes it possible to predict plant miRNA targets using a homology search. In this study, the predicted DE miRNAs were used to search for their potential targets against the cotton SGN Unigene Database (ftp://ftp.sgn.cornell.edu) using psRNATarget (http://bioinfo3.noble.org/psRNATarget) with default parameters. As a result, a total of 466 predicted target genes (Supplement Table 6) were identified that may be regulated by the 27 DE miRNAs. Among these predicted targets, some are known to be related to pollen and flower development and cell cycles. For example, SQUAMOSA Promoter-Binding Protein-Like (SPL-like genes), the predicted targets of miR156 which was down-regulated in CMS line, are important transcription factors involved in regulation of flowering and vegetative phase change [26,27,28,29,30]. TEOSINTE BRANCHED1, CYCLOIDEA, and PCF ( TCP4), the predicted target of up-regulated miR319 in the CMS-D2 line, is involved in stamen development [31]. Auxin response factor (ARF) genes, the predicted targets of up-regulated miR894 and PN ghr-60, and down-regulated miR167 in the CMS line, are essential for vegetative and reproductive development [32]. MYB genes, the predicted target genes of down-regulated miR159 and miR169 in the sterile genotype, are important transcription factors involved in pollen development [33]. Ribonucleotide reductase (RNR2A), the predicted target gene of up-regulated miR166 and miR165 in the CMS line, is critical for cell cycle progression and DNA damage repair [34]. D-type cyclins (CYCD2), the predicted target gene of up-regulated miR403 in the sterile line, is involved in G1 phase of mitotic cell cycle [35]. RETINOBLASTOMA RELATED (RBR), the predicted target gene of down-regulated miR167 in the CMS line, is involved in G1/S transition of mitotic cell cycle and regulation of DNA endoreduplication [36].

Based on the information of protein similarities supplied by NCBI (ftp://ftp.ncbi.nih.gov/repository/UniGene/), a Gene Ontology (GO) enrichment analysis was performed for the DE miRNA targets according to the Arabidopsis thaliana's GO information. For the 16 up-regulated miRNAs in the CMS line, their predicted target genes mainly participate in the process of unknown biological process, regulation of transcription, response to salt stress, metabolic process, embryo development ending in seed dormancy, and response to abscisic acid stimulus (Fig 2). Those targets of down regulated miRNAs mainly regulate the genes involved in unknown biological process, regulation of transcription, response to cadmiumion, oxidation-reduction process and transport, and metabolic process (Fig 3).

Quantitative RT-PCR

Since miR156 has been reported to target and regulate the SPL-like genes which is one of the important transcription factors related to pollen development, the expression profile of miR156 and its target SPL-like genes were chosen for further verification using quantitative RT-PCR (qRT-PCR). The result confirmed that miR156 was expressed at a lower level in the CMS line (Fig 4). Subsequently, using a conserved domain search in the cotton EST database of NCBI, 22 unique SPL-like genes (named as Ghspl1- Ghspl22) were obtained in cotton (Supplement sequences.doc). Based on sequence analysis using psRNATarget [37], 11 genes (Ghspl1, Ghspl3, Ghspl6, Ghspl7, Ghspl15, Ghspl16, Ghspl17, Ghspl18, Ghspl19, Ghspl20 and Ghspl21) with predicted miRNA/target pairs were considered as the target genes of miR156, while the other 11 genes were considered as nontarget genes. The qRT-PCR results revealed that 10 of the 22 SPL-like genes, including 4 target genes (Ghspl3, Ghspl7, Ghspl18, and Ghspl21) and 6 nontarget genes (Ghspl2, Ghspl4, Ghspl8, Ghspl11, Ghspl14, and Ghspl22) were differentially expressed between the CMS line and the fertile line. Seven of them were up-regulated (Ghspl2, Ghspl3, Ghspl4, Ghspl7, Ghspl8, Ghspl11 and Ghspl14) (Fig 5A) and 3 down-regulated (Ghspl18, Ghspl21 and Ghspl22) in the CMS line (Fig 5B).

Discussion

Of the annotated miRNAs from animals and plants, only a small portion have been identified in cotton by in-silico analysis [18,19,38] and deep sequencing [20], and little is known about their associations with pollen development. In this study, deep sequencing was performed to identify miRNAs from cotton flower buds of 3 mm in length (at roughly the stage of male meiosis) in CMS-D2 and its isogenic maintainer line. This resulted in the identification of 11 known cotton miRNAs, 32 conserved plant miRNAs and 86 predicted miRNAs. The results implied that the expression of miRNAs has tissue specificity and some of the novel miRNAs may be specifically expressed in cotton flowering buds.

According to the sequence reads, miR167, miR166 and miR172 were the most abundant miRNAs in the two genotypes. Previous studies [20,39] pointed out that some highly conserved miRNA families such as miR156/157, miR167 and miR172 have high levels of abundance which may represent a relationship between evolutionary conservation and expression abundance. Here, similar expression profile could be found for miR167 and miR172. For miR156, Kwak et al. [20] reported about one hundred thousand reads could be sequenced in a single library of ovules. In the current study, however, the expression abundance of miR156 was significantly lower (<1000 reads) in flowering buds of both the CMS and maintainer lines, especially in the CMS line. Thus, these results imply that some miRNAs are very important regulation factors during meiosis with tissue specific expressions, even though they were previously thought to have conserved expression abundance.

It is known that SPL-like genes participated in the biological processes of leaf development [40], cell number and size [41], plant architecture [42], GA signaling [43], response to copper and fungal toxin [44,45], flower and fruit development [26,27,28,29] and sporogenesis [30]. The current study further confirmed that SPL-like genes function in the early stage of anther and pollen development [46]. Yang et al. [47] reported that a probable pathway responsible for CMS involving mitochondrial retrograde regulation in stem mustard (Brassica juncea var. tumida), with SPL-like genes as a target nuclear gene for a mitochondrial signal. It was reported that miR156 mainly regulates the expression of some SPL-like genes by a post-transcriptional way [48]. In some studies, abnormal expression of miR156 and some SPL-like genes resulted in male sterility in plants [49,50]. This implied that the miR156-target and nontarget SPL-like transcription factors act in concert to ensure male fertility. In our study, the differential expression profiles of miR156, and its 4 target and 6 nontarget SPL-like genes were confirmed by quantitative RT-PCR between the CMS-D2 line and its isogenic maintainer line. In the 4 target SPL-like genes, only two (Ghspl3 and Ghsp7) were up-regulated in the CMS line, as expected. This discrepancy may be explained by the following reasons. miRNAs in plants tend to target the coding sequence (CDS), but some can be in the untranslated regions (UTR) [51,52]. With UTR modifications and splice variants, miRNA binding sites may change, leading to inactivation of the binding affinity for a given miRNA [53], as well as the overall number of miRNA binding sites for the gene [54]. Therefore, a down-regulated miRNA may not necessarily lead to up-regulation of its target mRNA. This may in part explain why some SPL-like genes remained unchanged or were even up-regulated.Furthermore, recent studies suggested that translation inhibition may be more common in plants than previously anticipated [55,56,57]. Therefore, just like AtSPL3 [58], Ghspl18 and Ghsp21 were possibly affected by translational inhibition in the CMS-D2 line. miR156-SPL-like genes could provide an interesting system for further studies in CMS cotton. As for those differentially expressed miR156-nontargert genes, except for Ghspl22 the other five genes were all up-regulated in the CMS line. Whether the over-expression of the notarget SPL-like genes in the CMS line is the complementary consequence of abnormal expression of some miR156 targeted SPL-like genes is unknown. It is likely that, as the target SPL-like genes for a mitochondrial signal, the miR156-SPL feedback interaction may be in part responsible for CMS-D2 in cotton involved in mitochondrial retrograde regulation.

Although almost all miRNAs were identified originally from sporophytes, most of them were also found to be expressed in developing pollen, including pollen-enriched ones [16,17,59,60,61,62]. Most miRNAs detected in Arabidopsis mature pollen were in a low abundance [15], possibly because the mature pollen is in the terminal development status ready for fertilization and needs no more miRNAs for gene regulation [17]. Therefore, the existence and high expression of some miRNAs in developing flowering buds (sporophytes) demonstrated their importance for male including pollen development. In this study, a total of 129 miRNAs were identified in the two genotypes, including some novel miRNAs not previously reported in cotton. Our results will further enrich the information of miRNAs related to male development in plants. By comparing the sequence reads of each miRNA identified in the two genotypes, a total of 27 (including 5 known cotton miRNAs, 17 conserved miRNAs and 5 novel miRNAs) differentially expressed miRNAs were identified between the CMS line and its isogenic maintainer line at the stage of male meiosis. Among the miRNAs identified, about half of (22/43) the previously reported miRNAs were differentially expressed in the two lines. For 86 novel miRNAs, only 5 were differentially expressed with low abundance. Therefore, the male sterile cytoplasm (CMS-D2) affected the expression profile of 27 miRNAs at the stage of male meiosis. As reported in Arabidopsis, there existed some pathways for different miRNAs [63]. For example, a miR156-AtSPL9-miR172 regulatory pathway was found to play a critical role for transition from juvenile to adult leaf development. Further studies are underway in our lab by using microspores and pollen at different stages of the CMS line and its isogenic maintainer line to detect the expression relationships between different miRNAs and their predicted target genes.

In the present study, 466 predicted targets were predicted for the 27 differentially expressed miRNAs, with an average of 17 target genes per miRNA. As mentioned by Pang et al. [21], the number of predicted targets per miRNA in cotton was larger than that in Arabidopsis, suggesting additional paralogous and homoeologous genes in the allotetraploid cotton species. Recent studies showed that precise expression of a complex set of transcription factors is crucial for plant development, especially for pollen development [64,65]. Many predicted targets of miRNAs encode for transcription factors, implying that some miRNAs are master regulators in plants [17]. In this research, based on the biological processes of the predicted target genes, the regulation of transcription accounted for the highest portion in known biological processes for both up- and down- regulated miRNAs in the CMS line. This further confirmed the critical role of miRNAs in male meiosis which leads to pollen development. As for the molecular function and cellular component, the differentially expressed miRNAs and their predicted target genes may be mostly involved in binding and function in cytosol and plasma membrane. These results indicate that the differentially expressed miRNAs affected by the male sterile cytoplasm may play a crucial role in nucleo-cytoplasmic interactions.

In addition to the SPL-like genes, some previously reported genes in other plants, including MYB, ARF, F-box and AP2 which are related to pollen development and flower development, were also identified as the predicted target genes of the differentially expressed miRNAs. Some reports have demonstrated that overexpression of miR156 [49], miR159 [66,67] and miR167 [68] repressed the mRNAs levels of some SPL-like genes, MYB genes and ARF genes, respectively, which affected pollen development, leading to male sterility. Furthermore, some predicted target genes (RBR1, CYCD2 and RNR2A) related to meiosis were also identified in the CMS line. RBR is required for correct differentiation of male gametophytic cell types at early meiotic events, and loss of RBR during meiosis resulted in a reduction of crossover formation and an associated failure in chromosome synapsis [69], leading to prevent or delay in cell determination during plant male gametogenesis [70]. Interestingly, RBR protein interacts specifically with D-type cyclins (CYCD) and regulates cell proliferation by binding and inhibiting E2F transcription factors [71]. These two genes were both identified as the predicted target genes of down-regulated miR167 and up-regulated miR403 respectively. RNR catalyzes a rate-limiting step in the production of dNTPs (deoxyribonucleotides) needed for DNA replication and repairing. Defective RNR often leads to cell cycle arrest, growth retardation, and p53-dependent apoptosis [72]. Mutants of RNR2 were more sensitive to UV-C light, and exhibited increased DNA damage, massive programmed cell death, and release of transcriptional gene silencing [73]. The sterility of CMS-D2 line was found to be sensitive to high temperature and light. Whether there exist some connections between this phenomenon and RNR2 needs further studies.

According to our previous report [6], the mitochondria atpA gene may be the possible gene conferring male sterility in the CMS-D2 system of cotton. In this study, some differentially expressed miRNAs (miR159, miR160, miR169, miR172, miR482, miR1311, and PN-ghr-60) that may participate in ATP pathway were also identified. Their predicted target genes are mainly involved in ATP binding, ATP transmembrane transporter activity, ATP:ADP antiporter activity, ATP-dependent DNA helicase activity at the mitochondrial inner or outer membrane. As noted previously, ATP production is very important for supplying of energy for normal meiosis and pollen development and and deficiency in ATP production was found in many CMS systems [74,75]. The results in the present study using CMS-D2 cotton provide evidence that genes involved in ATP production in CMS plants may be suppressed by miRNAs.

Conclusions

Using deep sequencing, 129 miRNAs were identified in the developing flowering buds of CMS-D2 line and its isogenic maintainer line, 27 of which were differentially expressed between the two lines. 466 putative targeted genes were predicted including transcription factors such as SPL-like genes, MYB, ARF, F-box and AP2, and genes involved in ATP production. The results indicate that the male sterile cytoplasm affect gene expression of miRNAs which in turn may down regulate their target genes in mitochondria-nuclear interactions. Further detailed studies using various male tissues and development stages will provide a better understanding of the molecular mechanisms of male sterility and mitochondria-nuclear interactions in Upland cotton with the CMS-D2 cytoplasm.

Methods

Plant materials

CMS cotton (CMS-D2) ZBA with cytoplasm from G. harknessii (D2 genome) and its isogenic maintainer ZB with Upland cotton (AD1) cytoplasm were provided by Cotton Research Institute (CRI), Chinese Academy of Agricultural Science (CAAS), Anyang, Henan, China. Plants were grown in the Cotton Research Farm, CRI-CAAS, near Anyang. According to previous results of cytological investigation for CMS cotton [76], flower buds about 3 mm in length (at roughly the stage of male meiosis) were collected from the field and used in present research. All harvested samples were snap-frozen in liquid nitrogen and stored at -80°C before use.

Small RNA Library Preparation and Sequencing

Total RNA was extracted using mirVana miRNA Isolation kit (Ambion, USA) according to the manufacturer's instructions. Deep sequencing was performed by LC Sciences (Houston, USA). In brief, 20 μg of total RNA were used for sequencing library preparation. The small RNA fractions with the length of 15-40 nt were isolated by polyacrylamide gel electrophoresis and with proprietary adaptors (Illumina, USA). The short RNAs were then converted to cDNA by RT-PCR [77]. The PCR products were size-fractionated and recovered for sequencing using Illumina's Genome Analyzer II according to manufacturer's instructions (Illumina, USA).

Data processing

The sequencing data were processed with an commercially available program called ACGT101-miR3.5 (LC Sciences, USA). The program parameters were set mainly according to the previous published paper [77] with some modifications for better miRNA prediction in plants. The sequenced sequences (sequ-seqs) that passed through a series of additional filters are called "mappable sequences", which are sequences satisfying the following criteria from the statistics of plant miRNAs in miRBbase 16.0: (1) not sequencing adapters; (2) containing no more than 80% A, C, G, or T; (3) containing no more than three N (undermined bases); (4) containing no stretches of A7, C8, G6, or T7; (5) observed more than two times; (6) being longer than 15 nt and shorter than 26 nt; and (7) not originated from known classes of RNAs (i.e., mRNA, rRNA, tRNA, snRNA, snoRNA and repetitive sequence elements).

The "mappable sequences" were compared with those miRNA precursor/mature miRNAs in the miRbase (http://www.mirbase.org/index.shtml, released 15.0) to identify conserved plant mature miRNAs in cotton and other plants, the EST and GSS sequences of cotton obtained from NCBI database (http://www.ncbi.nlm.nih.gov/) were used for mapping miRNAs and pre-miRNAs. For identification of potentially novel miRNAs, we followed the community standard of Meyers's report [78]. And RNA hairpin structure was predicted using RNAfold web server (http://rna.tbi.univie.ac.at/cgi-bin/RNAfold.cgi). A sequence to be considered a potential miRNA needs to form a perfect structure with parameters below, (1) Number of allowed errors in one bulge in stem <= 12; (2) Number of basepair (bp) in stem region >= 20; (3) Free energy (dG in kCal/mol) <=-17; (4) Length of hairpin (up and down stem + terminal loop) >= 60; (5) Length of terminal loop <= 350; (6) Number of allowed errors in one bulge in mature region <= 8; (7) Number of allowed biased errors in one bulge in mature region <= 2; (8) Number of allowed biased bulges in mature region <= 2; (9) Number of basepair (bp) in mature or mature region >= 15; (10) Percentage of small RNA in stem region (pm) >= 100%; and (11) Number of allowed errors in mature region <=6.

In-silico analysis of SPL-like genes

Previous study pointed out that sequence CX4CX13HX5HX15CQQCX3HX11C is a conserved domain of SPL-like gene [79]. Here this domain was used to conduct blastx search in NCBI (http://www.ncbi.nlm.nih.gov/) in the cotton EST database using the default parameters. Subsequently, primers were designed using software of Oligo 6 for each uniEST.

Quantitative RT-PCR

Similar to other reports on miRNAs in plants [15,80,81,82,83], real-time PCR reactions of cDNA were performed using SYBR® Premix Ex TaqTM (Perfect Real Time) (TaKaRa, Japan) according to the manufacturer's instructions. Briefly, 2 l of cDNA template and 800 nM of each primer were added to 10 l of the 2Ã-SYBR® Premix Ex TaqTM, then ddH2O was added to a final volume of 20 l. PCR analysis was performed using CFX96TM (Bio-Rad, U.S.A). The reaction was initiated with pre-denaturation at 95ï‚°C for 30 s, followed by 40 cycles of denaturation at 95ï‚°C for 5 s, annealing at 58ï‚°C for 20 s and extension at 72ï‚°C for 30 s. Immediately after the final PCR cycle, a melting curve analysis was done to determine specificity of the reaction. Reactions were typically done in triplicate and the controls without template were included for each gene. 5S rRNA and 18S rRNA were used as internal control for miRNA156 and SPL-like genes, respectively. The primers used in real-time PCR are shown in Supplement Table 7.

The threshold cycle (Ct) values of each reaction were determined automatically by the instrument software, and the relative amount of each gene to internal control was calculated using the equation 2-ΔΔCt, whereΔΔCt = (Ct target - Ct 5s/18s rRNA)sample X - (Ct target - Ct 5s/18s rRNA) sample1. The whole assay protocol was repeated three times to ensure the reliability of the assay data. The standard deviations of the data were obtained from the three independent experiments. Statistical significance of expression differences was analyzed using the Student t-test.