Transcriptome/Degradome-Wide Discovery of MicroRNAs and Transcript Targets in Two Paulownia australis Genotypes

Abstract

MicroRNAs (miRNAs) are involved in plant growth, development, and response to biotic and abiotic stresses. Most of the miRNAs that have been identified in model plants are well characterized, but till now, no reports have previously been published concerning miRNAs in Paulownia australis. In order to investigate miRNA-guided transcript target regulation in P. australis, small RNA libraries from two P. australis (diploids, PA2; and autotetraploids, PA4) genotypes were subjected to Solexa sequencing. As a result, 10,691,271 (PA2) and 10,712,733 (PA4) clean reads were obtained, and 45 conserved miRNAs belonging to 15 families, and 31 potential novel miRNAs candidates were identified. Compared with their expression levels in the PA2 plants, 26 miRNAs were up-regulated and 15 miRNAs were down-regulated in the PA4 plants. The relative expressions of 12 miRNAs were validated by quantitative real-time polymerase chain reaction. Using the degradome approach, 53 transcript targets were identified and annotated based on Gene Ontology and Kyoto Encyclopedia of Genes and Genomes pathway analysis. These targets were associated with development, stimulus response, metabolism, signaling transduction and biological regulation. Among them, 11 targets, including TCP transcription factors, auxin response factors, squamosa promoter-binding-like proteins, scarecrow-like proteins, L-type lectin-domain containing receptor kinases and zinc finger CCCH domain-containing protein, cleaved by four known miRNA family and two potentially novel miRNAs were found to be involved in regulating plant development, biotic and abiotic stresses. The findings will be helpful to facilitate studies on the functions of miRNAs and their transcript targets in Paulownia.

About the Authors:

Suyan Niu

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Guoqiang Fan

* E-mail: [email protected]

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Enkai Xu

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Minjie Deng

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Zhenli Zhao

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Yanpeng Dong

Affiliations Institute of Paulownia, Henan Agricultural University, Zhengzhou, Henan, P.R. China, College of Forestry, Henan Agricultural University, Zhengzhou, Henan, P.R. China

Introduction

MicroRNAs (miRNAs), a class of endogenous 21–24 nucleotide (nt), single-stranded, non-coding RNAs, modulate the gene expression at the transcription and post-transcriptional levels [1], [2]. In plants, emerging data demonstrate that miRNAs play vital regulatory roles in a wide range of biological processes, including regulation of plant growth, development, signal transduction [3]–[6], and response to biotic and abiotic stresses via interactions with their specific target mRNAs [7]–[10]. Thereby, understanding of miRNAs functions in plant requires the recognition of their target genes. Initially, based on the either perfect or near-perfect sequence complementarity between miRNAs and their target mRNAs, the plant miRNA transcript targets have been identified via computational prediction [11], [12]. However, computational target prediction method is very challenging because of the existence of a higher mismatch in miRNA-target mRNA pairing [13]. With the development of high-throughput sequencing technology, recently, a transcriptome wide experimental method combining high-throughput miRNA profiling with degradome sequencing analysis, has been successfully developed to identify miRNA-directed mRNA cleavage at a large scale [14], [15]. Using this strategy, so many plant species have been studied, for instance, Arabidopsis [15], Oryza sativa [16], maize [17], Brassica napus [18], Medicago truncatula [19], Fragaria ananassa [20], grapevine [21], Taxus [22], Populus trichocarpa [23], Populus euphratica [24], Trifoliate orange [25], Paulownia tomentosa [26], and Paulownia fortunei [27]. These discoveries have triggered detailed biological studies on gene regulation with the result that the number of reported plant miRNAs and their transcript targets has been increasing rapidly.

Paulownia is a genus of deciduous tree species in the family Paulowniaceae, related to and sometimes included in the family Scrophulariaceae. It is a very adaptable, extremely fast-growing woody plant, and in managed plantations, it can be harvested for saw timber in as little as 5 years. Paulownia can not only provide wood for a variety of products including timber, fuelwood, herbal medicines, boxes, clogs, musical instruments, and surfboards, but can also benefit the environment and ecology; for example, it can increase crop production when used for intercropping on farmland, and prevent soil erosion [28], [29]. Because of its specific characteristics and economic value, in its native China, Paulownia is popular for reforestation, intercropping on farmland, and roadside and garden planting. To enlarge the germplasm and breed new varieties, in recent years, several autotetraploid paulownia seedlings were induced successfully using colchicines [30]–[33]. The autotetraploid cultivars of paulownia contains two sets of the same chromosomes, and compared with their diploid progenitors, the autotetraploids have better timber quality and improved stress resistance [34], [35], which makes understanding the molecular mechanisms underlying the differences of diploids and autotetraploids imperative. To our knowledge, the different paulownia species have genetic diversity and differentiation, although miRNAs in the P. tomentosa [26] and P. fortunei species have been identified [27], there have been no detailed reports of miRNAs in P. australis. In the current study, we used Solexa sequencing and degradome technology to analyze four sequencing libraries constructed from the seedlings of the two P. australis genotypes to identify conserved and novel miRNAs and to investigate the potential roles of their transcript targets.

Materials and Methods

Plant culture and treatment

Plantlets of two P. australis genotypes (diploids, PA2; and autotetraploids, PA4) were cultured in vitro on medium as described previously [30], [36]. All the samples were grown at 25°C under a 16/8 h (light/dark) photoperiod for 30 days, and at least three parallel samples were prepared for each genotype. After culturing, samples were collected from nine different seedlings of each genotype and subsequently mixed, frozen immediately in liquid nitrogen, and stored at −86°C until further use.

Construction and sequencing of the small RNA library

Total RNA was isolated from the leaves of the PA2 and PA4 samples using TRIzol reagent (Invitrogen, Carlsbad, CA, USA) following the manufacturer's instructions. Two small RNA (sRNA) libraries were constructed, one from the PA2 sample and the other from the PA4 sample, using methods described elsewhere [24], [37]. Briefly, bands in the 18–30 nt size range were purified by electrophoretic separation on a 15% denaturing polyacrylamide gel and then successively ligated with 5′ and 3′ adapters. After reverse transcription and PCR, the amplified products were subjected to Solexa sequencing (Illumina, San Diego, CA, USA) at the Beijing Genomics Institute, Shenzhen, China. The data used in this publication have been deposited in the NIH Short Read Archive database (http://www.ncbi.nlm.nih.gov/sra) with accession number SRP041440 and SRP041442 (Alias: PRJNA245384 and PRJNA245383).

Identification of conserved and novel miRNAs in the P. australis libraries

Low-quality reads, poly(A) reads, oversized insertions, reads shorter than 18 nt, and adaptor contaminated reads were filtered out from the raw reads to yield the clean reads. The length distribution of the 18–30 nt clean reads was analyzed and the reads were mapped to the sequences in the P. australis UniGene database (http://www.ncbi.nlm.nih.gov/sra; Accession Number SRP032321) using SOAP [38]. Perfectly matched reads were analyzed by running Blastall (http://www.ncbi.nlm.nih.gov/staff/tao/URLAPI/blastall/) against the Rfam (http://www.sanger.ac.uk/software/Rfam) and GenBank databases (http://www.ncbi.nlm.nih.gov/) to discard tRNAs, rRNAs, snRNAs, snoRNAs, and other ncoRNAs. The remaining reads were annotated by aligning them against the sequences in the miRBase database 19.0 (http://www.mirbase.org/) to identify known miRNAs. When a read in the P. australis libraries shared homology with fewer than two mismatches with a miRNA sequence in miRBase, it was considered to be an evolutionarily conserved miRNA.

The MIREAP software (https://sourceforge.net/projects/mireap/) was used to detect novel miRNAs in the two libraries by predicting the stem-loop structures and estimating the minimal folding free energy of the unannotated sRNAs (no matches in miRBase) that mapped to the P. australis UniGene sequences. The stem-loop structures of the candidate precursor miRNAs (pre-miRNAs) were predicted using Mfold (http://mfold.rna.albany.edu/?q0mfold) [39]. The strict criteria used to annotate the candidate miRNAs were as described by Meyers et al. [40].

Differential expression analysis of miRNAs in the two P. australis libraries

To identify the differentially expressed miRNAs between the PA2 and PA4 libraries, the abundance of the miRNAs in the two libraries were normalized to one million, regardless of the total number of miRNAs in each sample. Then, statistical analysis was performed based on a Poisson distribution. Finally, when the fold-change was greater than 1.0 or less than −1.0 and the P-values were less than 0.05, the miRNA was considered to be significantly differentially expressed. The fold-changes and P-values [41], [42] were calculated from the normalized expression using the following formulas:

Normalized expression  =  (actual miRNA reads/total count of clean reads) × 1,000,000

Fold-change  =  log2 (PA4 normalized reads/PA2 normalized reads)

P-values were calculated as,

Degradome library construction, data analysis, and target identification

Two degradome libraries were constructed from the leaves of the PA2 and PA4 plants based on a protocol published previously [14], [15]. Briefly, using T4 RNA ligase (Takara, Dalian, China), poly(A) enriched RNA was ligated to a 5′RNA adapter containing a MmeI recognition site. Reverse transcription was performed to generate first-strand cDNA, followed by PCR enrichment, and digestion using the MmeI restriction enzyme (NEB, Ipswich, MA, USA). A double-stranded DNA adapter was then ligated to the digested products using T4 DNA ligase (NEB, Ipswich, MA, USA) and gel purified for PCR amplification. The resulting PCR products were sequenced using an Illumina HiSeqTM 2000 system. Low-quality sequences and adapters were removed from the raw reads of the degradome libraries, then the remaining reads were aligned to the P. australis UniGene database (http://www.ncbi.nlm.nih.gov/sra; Accession Number SRP032321) using SOAP software (http://soap.genomics. org.cn/) to define the coverage rate. The PAIRFINDER software (version 2.0) [22] was used to identify the miRNA-mediated cleaved fragments. Alignments with scores not exceeding four and no mismatches at the cleavage site (between the 10th and 11th nucleotides) were considered to be the potential miRNA targets. Furthermore, t-plots were built according to the distribution of signatures (and abundances) along these transcripts to help analyze the miRNA targets and RNA degradation patterns. Based on the locations of the target genes in the P. australis UniGene sequences, putative target genes were selected manually and subsequently mapped to the previously identified genes that were annotated according to the annotations of their homologous sequences in the GenBank Nr and Nt databases, and in the Swiss-Prot database using BLASTX searches with an E-value cutoff of 10⁻⁵.

Verification of miRNAs and their targets by quantitative real-time polymerase chain reaction

Quantitative real-time polymerase chain reaction (qRT-PCR) was used to validate and measure the expressions of selected miRNAs and their targets obtained from the Solexa sequencing and degradome analysis. The qRT-PCR was performed according to the protocol described previously [43]. For the qRT-PCR, total RNA were isolated from the leaves of two individual PA2 and PA4 plants at two developmental stages (30-day-old in vitro plantlets, and two-year-old saplings from field). Three biological replicates of each stage were used. The qRT-PCRs were performed in triplicate on a CFX96TM Real-Time PCR System (Bio-Rad, Hercules, CA, USA) using a SuperScript III platinum SYBR Green one-step qRT-PCR kit (Invitrogen, Carlsbad, CA, USA) according to the manufacturer's instructions. The PCR conditions were 50°C for 3 min, 95°C for 5 min, then 40 cycles of 95°C for 15 s, 55°C for 30 s, and 40°C for 10 min. Specific stem-loop primers and other primers for the miRNAs were designed based on the mature miRNA sequence. The U6 was used as an endogenous control. All the primers used for the qRT-PCR are listed in Table S1. The primers for the target genes were designed using Beacon Designer, version 7.7 (Premier Biosoft International, Ltd., Palo Alto, CA, USA), and the 18S rRNA of Paulownia was chosen as an endogenous reference gene for normalization. The relative expression level of a miRNA or target gene was calculated according to the method described previously [44].

Results

Analysis of sRNAs in the two P. australis libraries

A total of 13,895,340 and 13,537,466 raw reads were obtained in the PA2 and PA4 sRNA libraries, respectively. After removing low-quality reads, poly(A) reads, oversized insertions, reads shorter than 18-nt, and adaptor contaminated reads, we obtained 10,691,271 (PA2) and 10,712,733 (PA4) clean reads that represented 2,006,153 (PA2) and 2,418,971 (PA4) unique reads with lengths ranging from 18 to 30-nt. The majority of clean reads, 63.24% (PA2) and 68.79% (PA4), ranged from 20 to 24-nt in length, among which the 24-nt long reads were the most abundant, followed by the 21-nt long reads (Figure S1). The abundances of the 24-nt and 21-nt long sRNAs in the PA4 library were 6.48% and 0.04% more, respectively, than their abundances in the PA2 library. This atypical situation has also been reported in other hardwood species, such as Liriodendron chinense [45].

The clean reads were used to query the P. australis UniGene database (http://www.ncbi.nlm.nih.gov/sra), the non-coding RNAs sequences deposited in GenBank (http://www.ncbi.nih.gov/Genbank), the Rfam database (http://rfam.sanger.ac.uk/), and miRBase 19.0 (http://microrna.sanger.ac.uk/sequences). These searches allowed us to assign annotations to each sRNA sequence, including rRNA, tRNA, snRNA, snoRNA, and known miRNA (Table S2). As a result, 848,652 (PA2) and 931,989 (PA4) known miRNA reads were detected. Finally, 1,809,478 (PA2) and 2,308,527 (PA4) unique unannotated sRNAs remained for predicting potentially novel miRNAs.

Identification of conserved miRNAs in P. australis

To detect conserved miRNAs in the two libraries, the unique reads were compared with the known miRNAs in miRBase 19.0, allowing two mismatches. In the two libraries, a total of 45 conserved miRNAs belonging to 15 miRNA families were found to share high identity with known plant miRNAs (Table 1 and Table S3). Among these 45 miRNAs, 40 and 43 were identified in PA2 and PA4, respectively; 38 of them were expressed in both libraries, two were found only in the PA2 library, and five occurred only in the PA4 library. Most of the conserved miRNAs (91.67%) were 21-nt long and the remaining miRNAs were 20 or 23-nt long (Table 1 and Table S3). This result is consistent with the current understanding that canonical miRNAs are 21-nt long, while canonical small interfering RNAs (siRNAs) are 24-nt long [46]. Among the 15 miRNA families, the miR166 family was dominant in both libraries, accounting for 89.44% (PA2) and 86.35% (PA4) of all conserved miRNA reads, followed by the miR159 family (Table 1). Several miRNA families, such as miR156, miR396, miR397, miR398, and miR482, had moderate expression levels in the two libraries, while other miRNA families showed very low levels of expression, with fewer than 100 reads, in both libraries. Moreover, different members in the same miRNA family displayed significantly different expression levels. For instance, members of the miR166 family varied in abundance from 50,932 to 563,946 reads in the two libraries (Table 1).

[Figure omitted. See PDF.]

Table 1. Conserved miRNAs identified from P. australis.

https://doi.org/10.1371/journal.pone.0106736.t001

Identification of novel miRNAs in the P. australis sRNA libraries

The characteristic stem-loop structure of pre-miRNA was employed to predict novel miRNAs [40]. We found that 31 of the potential pre-miRNAs met this requirement; 15 pre-miRNAs exhibited the 3p:5p miRNA pair, providing more evidence that they were novel miRNAs, while the other 16 miRNAs were also considered to be potentially novel miRNAs (Table 2 and Table S4). Among the 31 candidate novel miRNAs, 10 miRNAs were common between the two libraries, while 10 and 11 were specific to PA2 and PA4, respectively. The length of the mature miRNAs varied from 20 to 23 nt, with the majority being 21-nt long. The mature miRNA sequences were localized inside the stem-loop structures, with almost half in either the 3p or 5p arms. We observed that the average pre-miRNA length was 151 nt, and these precursors had minimal folding free energies that ranged from −142.5 kcal mol⁻¹ to −30.8 kcal mol⁻¹ with an average of −56.3 kcal mol⁻¹. In this study, compared with the abundance of the conserved miRNAs, the majority of novel miRNAs had relatively low expressions, and only five (PA2) and six (PA4) of the novel miRNA candidates had more than 1,000 reads (Table 2 and Table S4).

[Figure omitted. See PDF.]

Table 2. Novel miRNAs identified from P. australis.

https://doi.org/10.1371/journal.pone.0106736.t002

Expression patterns of conserved and novel miRNAs

The differential expression analysis was performed based on the normalized read counts for each identified miRNA. A total of 41 miRNAs showed statistically significant changes (fold change ≥1.0 or ≤−1.0, and P-values ≤0.05) in their relative abundance between the two libraries. Among them, 26 miRNAs were up-regulated (13 conserved miRNAs and 13 novel miRNAs; fold-change ≥1.0) and 15 were down-regulated (three conserved miRNAs and 12 novel miRNAs; fold-change ≤−1.0) in the PA4 library compared with in the PA2 library. The expression levels of some miRNAs changed significantly. For example, the expression levels of pas-miR169b-3p, pas-miR169c-3p, pas-miR396c-3p, pas-miR396d-3p, pas-miR171a, pas-miR171b, and pas-miR171c increased or decreased by about 5-fold in the PA4 library (Table 1). We also detected 21 novel miRNAs in the PA4 library that had 5-fold greater or lesser expression levels than in the PA2 library (Table 2).

Identification of miRNAs transcript targets in P. australis by using degradome analysis

To better understand the physiological functions and biology processes that these miRNAs may be involved in during the development of P. australis, transcript target identification was performed based on the degradome approach. Using PAIRFINDER, a total of 53 targets were predicted to be cleaved by 11 of the conserved miRNA families and three novel miRNA candidates (Table S5). The target transcripts were pooled and grouped into three categories (I–III) according to their relative abundances [14], [15] (Figure 1). Among these identified targets, 32 (37 cleavage sites) belonged to category I, and 20 (34 cleavage sites) and six (six cleavage sites) targets belonged to categories II and III, respectively (Table S5). We then performed a Gene Ontology (GO) analysis to assign functional annotations to the predicted target genes, as described previously [47]. More than 80% of the target genes were annotated as being involved in regulation of biological processes, and the GO annotations for the predicted target genes are shown in Figure 2. The Kyoto Encyclopedia of Genes and Genomes (KEGG) Pathway database was used to further classify the miRNA target genes. The KEGG Pathway annotation showed that the target genes were involved in metabolic pathways, plant-pathogen interaction, pyruvate metabolism, carbon fixation in photosynthetic organisms, purine metabolism, RNA polymerase, pyrimidine metabolism, nitrogen metabolism, plant hormone signal transduction, cellular metabolism, and disease. (Table S5). The prediction and annotation of the miRNA target genes may provide some new insights into how P. australis miRNAs regulate gene expression in this plant.

[Figure omitted. See PDF.]

Figure 1. Target plots (t-plots) of miRNA targets in different categories confirmed by degradome sequencing.

(A) T-plot (top) and miRNA: mRNA alignments (bottom) for two category I targets, Unigene11776_All and CL1057.Contig2_All transcripts. The arrow indicates signatures consistent with miRNA-directed cleavage. The solid lines and dot in miRNA: mRNA alignments indicate matched RNA base pairs and GU mismatch, respectively, and the red letter indicates the cleavage site. (B) CL5129.Contig3_All and CL6407.Contig9_All, a category II target for pas-miR156d and pas-miR167. (C) CL11907.Contig1_All and Unigene1681_All, a category III target for pas-miR168b and pas-miR169e.

https://doi.org/10.1371/journal.pone.0106736.g001

[Figure omitted. See PDF.]

Figure 2. Gene Ontology analysis of miRNA targets in P. australis.

https://doi.org/10.1371/journal.pone.0106736.g002

Validation of the candidate miRNAs and their targets by qRT-PCR

To verify the existence and expression levels of the miRNAs determined from the high-throughput Solexa sequencing, 12 miRNAs with different expression levels were randomly selected for qRT-PCR analysis. As shown in Figure S2, the expression patterns of the miRNAs obtained by qRT-PCR showed similar trends to their expression patterns in the two libraries as determined by Solexa sequencing. Compared with their expressions in PA2, the expression levels of pas-miR156c, pas-miR398a-3p, pas-miR408a-5p, and pas-mir22-3p were up-regulated in PA4 in the 30-day plantlets and in the 2-year saplings, pas-miR319a-3p were on the contrary. Further, the expressions of pas-miR160a, pas-miR167, pas-miR171a, pas-miR397a, pas-mir1, and pas-mir14 in PA4 were up-regulated in the 30-day plantlets, and down-regulated in the PA4 in the 2-year saplings, while the expressions of pas-mir3 were the opposite to this. Thus, with the development of the plants, the expression levels of some of the miRNAs showed different trends. Six miRNAs (pas-miR160a, pas-miR167, pas-miR319a-3p, pas-miR398a,-3p pas-miR408a-3p, pas-mir1) had the same expression trend between PA2 and PA4 in the two stages. Furthermore, to confirm the reliability of degradome sequencing technology and the potential correlation between miRNAs and their transcript targets, 12 genes from two P. australis genotypes were also selected for qRT-PCR assays. The results showed that except the targets of auxin response factor (ARF) ARF8 (CL4211.Contig3) and scarecrow-like protein (SCL) SCL15 (CL10503.Contig1), the expression levels of the rest targets were inversely correlated with these of the corresponding miRNAs (Figure 3). During the different developmental stages, pas-miR160a, pas-miR167, pas-miR171 and pas-mir1 at a relatively higher level in the PA4 than in the PA2 at the 30-day plantlets stage, and lower level at the 2-year saplings stage, while its transcript targets, CL3173.Contig7, CL11603.Contig1, CL6407.Contig9, CL11078.Contig2, CL11078.Contig3, and Unigene9061 expressed in the reverse way as expected, and these coding for proteins are members of the auxin response factors ARF10, ARF18 and ARF6, the scarecrow-like proteins SCL6 and SCL22, and serine/threonine protein kinase, respectively. The expression levels of pas-miR319a-3p in PA4 were significantly lower than in the PA2 at two treatment stages, while the reverse was true for the transcription factor TCP4 (CL9103.Contig3) (Figure 3). Moreover, a reverse trend was noted between pas-miR156c and its target genes coding for squamosa promoter-binding-like protein (SPL) SPL6 and SPL12 (CL11428.Contig2 and CL5129.Contig2), and between pas-mir22-3p and its target gene coding for zinc finger CCCH domain-containing protein 53 (CL1197.Contig2) in the PA4 as compared to the PA2 (Figure 3). These results indicated that the miRNA and transcript target expression patterns were very complex and varied during the growth and development of P. australis. Furthermore, the possible roles of these miRNAs in the genome duplication changes from diploid to autotetraploid in these plants were revealed.

[Figure omitted. See PDF.]

Figure 3. Relative expression levels of the target genes in P. australis.

PA2-1, 30-day-old diploid in vitro plantlets; PA4-1, 30-day-old autotetraploid in vitro plantlets; PA2-2, two-year-old diploid saplings; PA4-2, two-year-old autotetraploid saplings. CL11428.Contig2 (SPL6) and CL5129.Contig2 (SPL12) targeted by pas-miR156c, CL3173.Contig7 (ARF10) and CL11603.Contig1 (ARF18) targeted by pas-miR160a, CL4211.Contig3 (ARF8) and CL6407.Contig9 (ARF6) targeted by pas-miR167, CL11078.Contig2 (SCL6), CL10503.Contig1 (SCL15) and CL11078.Contig3 (SCL22) targeted by pas-miR171a, CL9103.Contig3 (Transcription factor TCP4) targeted by pas-miR319a-3p; Unigene9061 (Serine/threonine protein kinase) targeted by pas-mir1; CL1197.Contig2 (Zinc finger CCCH domain-containing protein 53) targeted by pas-mir22-3p. Three independent biological replicates were performed. Values are means ± SD (n = 3).The expression levels of targets were normalized to 18SrRNA. The normalized miRNA levels in the PA2-1 were arbitrarily set to 1. *: Statistically significant differences between PA2 and PA4 under the same developmental stages (p-value was less than 0.05)

https://doi.org/10.1371/journal.pone.0106736.g003

Discussion

Variations in plant morphology and physiology resulting from genome duplication have occurred in many plants, such as Triticum, Gossypium, Spartina, Tragopogon, Brassica and Solanum [48], [49]. Generally speaking, genome duplication has led to the production of fast-growing, high-quality plants [50], [51]. The autotetraploid P. australis, which contains two sets of the same chromosomes, displayed apparent alterations in morphology, growth development, physiology, and gene expression when compared with their diploid counterparts [34], [35]. However, the mechanisms for these changes are poorly understood. MiRNAs are a class of endogenous small RNAs that have been involved in many processes, including growth, development, and resistance to stress and disease, by their ability to regulate gene expression in plants [52]–[54]. To understand the functions of miRNAs in diploid and autotetraploid P. australis, in this present study, we used Solexa sequencing and degradome approaches to construct two sRNA libraries and two degradome libraries from the PA2 and PA4 plants to identify conserved and novel miRNAs and their transcript targets. A total of 45 conserved miRNAs belonging to 15 miRNA families and 31 potential novel miRNA candidates along with 53 transcript targets were identified across the PA2 and PA4 libraries. Most of the identified conserved miRNA families are also conserved in other plant species, including Populus tomentosa [53], Populus euphratica [55], Oryza sativa [3], and Arabidopsis thaliana [56]. In the PA2 and PA4 libraries, the expression patterns varied dramatically among the different miRNA families, and different members in the same miRNA family also displayed significantly different expression levels. For instance, the read number varied from 26 (miR171 family) to 1,422,532 (miR166 family) (Table 1), and members of the miR166 family varied in abundance from 50,932 to 563,946 reads (Table 1). Moreover, the majority of novel miRNAs had relatively low expressions, and only five (PA2) and six (PA4) of the novel miRNA candidates had more than 1,000 reads (Table 2). Our results are in accordance with previous reports that novel miRNAs were often represented in relatively lower levels than conserved miRNAs [57], [58]. Furthermore, it is possible that the low-expression novel miRNAs may play particular functions in specific tissues, during developmental stages, or under various growth conditions. Whether these low-expression miRNAs are expressed at higher levels in other tissues or developmental stages, or are regulated by environmental stresses, remain to be investigated.

Previous studies have shown that many of the genes appeared to be methylated in tetraploid Paulownia plants specifically after genome duplication by the methylation-sensitive amplified polymorphism analysis [59]. DNA methylation has been reported to be involved in inducing gene silencing, which can restart or change the genes expression levels [60]. Salmon et al. [61] found that significant changes in DNA methylation patterns could explain the morphological plasticity and larger ecological amplitude of Spartina allopolyploids. Indeed, we found that the expression levels of many of the differentially expressed miRNAs in the PA4 library were not increased by more than two-fold compared with their expressions in the PA2 library. However, the expressions of about half of these miRNAs were significantly different in the two libraries. Some of the miRNAs were expressed at similar level in the PA4 and PA2 libraries (Table 2). These findings are similar to those reported previously in P. tomentosa, P. fortunei and tetraploid Arabidopsis thaliana lines [26], [27], [62], suggesting that the genome merger in the PA4 plants lead to nonadditive expression of the miRNA primary transcripts and miRNA target genes. Furthermore, the expression patterns of the differentially expressed miRNAs and transcript targets at different development stages were validated by qRT-PCR. The result showed that the differentially expressed miRNAs caused different expression levels in their transcript targets. Interestingly, we also found the expression levels of two target genes (CL4211.Contig3 and CL10503.Contig1) were inconsistent with those of their corresponding miRNAs (pas-miR167 and pas-miR171a). The similar phenomena were also observed in P. tomentosa, Phalaenopsis Aphrodite and cotton [26], [63], [64], indicating that the other mechanisms of regulating expression of the target genes exist. Above all, these results imply that the miRNAs with significantly varied expressions in the PA2 and PA4 P. australis plants are probably involved in the epigenetic changes of PA4 plants; however, the relation between the miRNA expression patterns and genome duplication may be more complex than we first thought.

To understand the biological functions of miRNAs, it is necessary to identify their transcript targets. In the present study, to avoid false-positive predictions of miRNA transcript targets in P. australis, we identified the 53 transcript targets for 11 miRNA families and three novel miRNA candidates by degradome sequencing, which opens up a new avenue for high-throughput validation of splicing targets [14], [15]. The target genes predicted for the conserved miRNAs in P. australis were similar or functionally related to validated plant miRNA targets, which were annotated as being involved in diverse physiological processes. For instance, pas-miR156 targeted the SPL protein family, which can affect diverse developmental processes such as leaf development, shoot maturation, phase change, and flowering in plants [65]; and pas-miR167 targeted ARF6 and ARF 8, which belong to a class of transcription factors known to control multiple processes in plants, including the regulation of gynoecium and stamen maturation, and seed dispersal [66], [67]. Thus, our results support the idea that conserved miRNAs take part in essential physiological processes in plants.

The analysis of the target genes identified for the differentially expressed miRNAs revealed that some of the target genes may play important roles in plant morphology and physiology. We found that the expression level of pas-miR319a-3p decreased by about 3-fold in PA4 compared with its level in PA2, and pas-miR319a-3p was predicted to target TCP transcription factors, which are plant-specific transcription factors that have been shown to participate in specifying plant morphological traits, such as organ border delimitation [68], cell division and proliferation [69], flower and leaf shape, and shoot outgrowth [70]. In Arabidopsis, the TCP transcription factors have been related to control of the morphology of shoot lateral organs and formation of the shoot meristem-dependent regulation of the expression of boundary-specific genes [71]. The putative transcription factor SCL is a number of GRAS protein family that are involved in several aspects of plant growth and development, including control of asymmetric cell division, maintenance of stem cell status, and induction of the regeneration of the root tip after laser ablation, and SCL expression has been associated with auxin distribution in root apical meristem [72]–[74]. In this study, three SCL genes were predicted to be targeted by the miR171 family, which was up-regulated in the PA4 plants. This result suggested that SCL may be involved in the formation of adventitious root and the other aerial organs of Paulownia. L-type lectin-domain containing receptor kinases (lecRKs) and zinc finger CCCH domain-containing protein were predicted to be targeted by pas-miR1 and pas-miR22-3p, respectively, suggesting that, besides their possible involvement in plant development, some of the miRNA target genes identified in this study could play fundamental roles in biotic and abiotic stresses. A previous study showed that LecRKs were most likely plasma membrane proteins, and were probably involved in mediating protein-protein interactions with a wide range of functions such as recognition of oligosaccharide or lipochitooligosaccharide signals, linking ABA-signaling, and response to salt, drought, cold stress, wounding, and disease in plants [75], [76]. Zinc finger CCCH domain-containing protein is a kind of RNA-binding protein and many studies have shown that it may be regulated by abiotic or biotic stresses, and could have regulatory functions in mRNA processing [77]–[80], thus supporting the possible roles for pas-miR1 and pas-miR22-3p in the adaptive response of PA4 to abiotic stress.

Based on the GO analysis of the targets of the identified miRNAs in the two P. australis genotypes plants, the target genes were separated into the biological process, cellular component, and molecular function. Some of the genes were annotated as being involved in biological regulation, cellular process, developmental process, response to stimulus, metabolic process, cell, cell part, organelle components, binding, and catalytic activity (Figure 2). Thus, our results suggest that the transcript targets might be closely related to the observed differences of phenotype (environmental adaptations) and resistance to biotic and abiotic stresses between the PA4 and PA2 plants. The functional role of the differentially expressed miRNAs will be the focus of future investigations. No target genes were predicted for many of the miRNAs identified in the two P. australis sRNA libraries, and a few of the predicted target genes were annotated as being of unknown function and hypothetical genes. Careful analysis of these potential targets will contribute further to our understanding of the role of miRNAs in P. australis.

In summary, miRNA target genes were identified using a degradome approach that included functional annotation and pathway analyses based on the GO and KEGG databases. Some of the transcript targets regulated by the differentially expressed miRNAs were related to the physiology and environmental adaptations. Our results suggest that the significantly varied expression miRNAs in the two P. australis genotypes are probably involved in the epigenetic changes of PA4 plants, the correlation between the miRNA expression pattern and genome duplication may be more complex than we first thought. Taken together, this study provides beneficial information for elucidating the miRNA-mediated regulation of transcript targets in P. australis and related species.

Supporting Information

[Figure omitted. See PDF.]

Figure S1.

Length distribution of sRNAs in P. australis. (A) Size distribution of total sequences. (B) Size distribution of unique sequences.

https://doi.org/10.1371/journal.pone.0106736.s001

(TIF)

Figure S2.

Results from qRT-PCR of miRNAs in P. australis. PA2-1, 30-day-old diploid in vitro plantlets; PA4-1, 30-day-old autotetraploid in vitro plantlets; PA2-2, two-year-old diploid saplings; PA4-2, two-year-old autotetraploid saplings. Three independent biological replicates were performed. Values are means ± SD (n = 3). The expression levels of miRNAs were normalized to U6. The normalized miRNA levels in the PA2-1 were arbitrarily set to 1. *: Statistically significant differences between PA2 and PA4 under the same developmental stages (p-value was less than 0.05).

https://doi.org/10.1371/journal.pone.0106736.s002

(TIF)

Table S1.

Quantitative RT-PCR validated miRNAs and their targets primers.

https://doi.org/10.1371/journal.pone.0106736.s003

(XLS)

Table S2.

Categories and statistical summary of sRNAs in P. australis.

https://doi.org/10.1371/journal.pone.0106736.s004

(DOCX)

Table S3.

The sequences of conserved miRNAs identified from P. australis.

https://doi.org/10.1371/journal.pone.0106736.s005

(XLSX)

Table S4.

The sequences of novel miRNAs identified from P. australis.

https://doi.org/10.1371/journal.pone.0106736.s006

(XLSX)

Table S5.

Identified targets of miRNAs involved in P. australis by degradome analysis.

https://doi.org/10.1371/journal.pone.0106736.s007

(XLSX)

Author Contributions

Conceived and designed the experiments: GQF SYN EKX MJD ZLZ YPD. Performed the experiments: GQF SYN EKX MJD ZLZ YPD. Analyzed the data: GQF SYN EKX MJD ZLZ YPD. Contributed reagents/materials/analysis tools: GQF SYN EKX MJD ZLZ YPD. Contributed to the writing of the manuscript: SYN GQF.

Citation: Niu S, Fan G, Xu E, Deng M, Zhao Z, Dong Y (2014) Transcriptome/Degradome-Wide Discovery of MicroRNAs and Transcript Targets in Two Paulownia australis Genotypes. PLoS ONE 9(9): e106736. https://doi.org/10.1371/journal.pone.0106736

References

1. Vaucheret H (2006) Post-transcriptional small RNA pathways in plants: mechanisms and regulations. Genes Dev 20: 759–771.

2. Meng Y, Shao C, Wang H, Chen M (2011) The regulatory activities of plant microRNAs: a more dynamic perspective. Plant Physiol 157: 1583–1595.

3. Sunkar R, Zhou X, Zheng Y, Zhang W, Zhu JK (2008) Identification of novel and candidate miRNAs in rice by high throughput sequencing. BMC Plant Biol 8: 25.

4. Bartel DP (2004) MicroRNAs: genomics, biogenesis, mechanism, and function. Cell 116: 281–297.

5. Voinnet O (2009) Origin, biogenesis, and activity of plant microRNAs. Cell 136: 669–687.

6. Mallory A, Vaucheret H (2006) Functions of microRNAs and related small RNAs in plants. Nat Genet 38: 31–36.

7. Sunkar R, Chinnusamy V, Zhu J, Zhu JK (2007) Small RNAs as big players in plant abiotic stress responses and nutrient deprivation. Trends Plant Sci 12: 301–309.

8. Shukla LI, Chinnusamy V, Sunkar R (2008) The role of microRNAs and other endogenous small RNAs in plant stress responses. Biochim Biophys Acta 1779: 743–748.

9. Chuck G, Candela H, Hake S (2009) Big impacts by small RNAs in plant development. Curr Opin Plant Biol 12: 81–86.

10. Zhou Z, Huang S, Yang Z (2008) Bioinformatic identification and expression analysis of new microRNAs from Medicago truncatula. Biochem Biophys Res Commun 374: 538–542.

11. Wang X, Reyes J, Chua N, Gaasterland T (2004) Prediction and identification of Arabidopsis thaliana microRNAs and their mRNA targets. Genome Biol 5: R65.

12. Adai A, Johnson C, Mlotshwa S, Archer-Evans S, Manocha V, et al. (2005) Computational prediction of miRNAs in Arabidopsis thaliana. Genome Res 15: 78–91.

13. Jones-Rhoades M, Bartel D (2004) Computational identification of plant microRNAs and their targets, including a stress-induced miRNA. Mol Cell 14: 787–799.

14. Addo-Quaye C, Eshoo TW, Bartel DP, Axtell MJ (2008) Endogenous siRNA and miRNA targets identified by sequencing of the Arabidopsis degradome. Curr Biol 18: 758–762.

15. German MA, Pillay M, Jeong DH, Hetawal A, Luo S, et al. (2008) Global identification of microRNA-target RNA pairs by parallel analysis of RNA ends. Nat Biotechnol 26: 941–946.

16. Zhou M, Gu L, Li P, Song X, Wei L, et al. (2010) Degradome sequencing reveals endogenous small RNA targets in rice (Oryza sativa L. ssp. indica). Front Biol 5: 67–90.

17. Shen Y, Jiang Z, Lu S, Lin H, Gao S, et al. (2013) Combined small RNA and degradome sequencing reveals microRNA regulation during immature maize embryo dedifferentiation. Biochem Biophys Res Commun 441: 425–430.

18. Zhou Z, Song J, Yang Z (2012) Genome-wide identification of Brassica napus microRNAs and their targets in response to cadmium. J Exp Bot 63: 4597–4613

19. Zhou Z, Zeng H, Liu Z, Yang Z (2012) Genome-wide identification of Medicago truncatula microRNAs and their targets reveals their differential regulation by heavy metal. Plant Cell Environ 35: 86–99.

20. Xu X, Yin L, Ying Q, Song H, Xue D, et al. (2013) High-throughput sequencing and degradome analysis identify miRNAs and their targets involved in fruit senescence of Fragaria ananassa. PLoS One 8: e70959.

21. Pantaleo V, Szittya G, Moxon S, Miozzi L, Moulton V, et al. (2010) Identification of grapevine microRNAs and their targets using high-throughput sequencing and degradome analysis. Plant J 62: 960–976.

22. Hao DC, Yang L, Xiao PG, Liu M (2012) Identification of Taxus microRNAs and their targets with high-throughput sequencing and degradome analysis. Physiol Plant 146: 388–403.

23. Shuai P, Liang D, Zhang Z, Yin W, Xia X (2013) Identification of drought-responsive and novel Populus trichocarpa microRNAs by high-throughput sequencing and their targets using degradome analysis. BMC Genomics 14: 233.

24. Li B, Qin Y, Duan H, Yin W, Xia X (2011) Genome-wide characterization of new and drought stress responsive microRNAs in Populus euphratica. J Exp Bot 62: 3765–3779.

25. Zhang JZ, Ai XY, Guo WW, Peng SA, Deng XX, et al. (2012) Identification of miRNAs and their target genes using deep sequencing and degradome analysis in trifoliate orange [Poncirus trifoliata L. Raf]. Mol Biotechnol 51: 44–57.

26. Fan G, Zhai X, Niu S, Ren Y (2014) Dynamic expression of novel and conserved microRNAs and their targets in diploid and tetraploid of Paulownia tomentosa. Biochimie 102: 68–77.

27. Niu S, Fan G, Zhao Z, Deng M, Dong Y (2014) High-throughput sequencing and degradome analysis reveal microRNA differential expression profiles and their targets in Paulownia fortunei. Plant Cell Tiss Org DOI 10.1007/s11240-014-0546-9.

28. Chirko CP, Gold MA, Nguyen PV, Jiang J (1996) Influence of direction and distance from trees on wheat yield and photosynthetic photon flux density in a Paulownia and wheat intercropping system. Forest Ecol Manag 83: 171–180.

29. Bayliss K, Saqib M, Dell B, Jones M, Hardy GSJ (2005) First record of ‘Candidatus Phytoplasma australiense’ in Paulownia trees. Australas Plant Path 34: 123–124.

30. Fan G, Wei Z, Yang Z (2009) Induction of autotetraploid of Paulownia australis and its in vitro plantlet regeneration. J Northwest Sci-Tech Univ 37: 83–90.

31. Fan G, Cao Y, Zhao Z, Yang Z (2007) Induction of Autotetraploid of Paulownia fortunei. Sci Silv Sin 43: 31–35.

32. Fan G, Yang Z, Cao Y, Q Z (2007) Induction of Autotetraploid of Paulownia tomentosa. Plant physiol Commun 43: 109–111.

33. Fan G, Zhai X, Wei Z, Yang Z (2010) Induction of autotetrapioid from the somatic cell of Paulownia tomentosa×Paulownia fortunei and its in vitro plantley regeneration. J Northeast Forestry Univ 38: 22–26.

34. Zhang X, Zhai X, Fan G, Deng M, Zhao Z (2012) Observation on Microstructure of Leaves and Stress Tolerance Analysis of Different Tetraploid Paulownia. J Henan Agric Univ 46: 646–650.

35. Zhai X, Zhang X, Zhao Z, Deng M, Fan G (2012) Study on Wood Physical Properties of Tetraploid Paulownia fortunei. J Henan Agric Univ 46: 651–654.

36. Fan G, Zhai X, Liu X (2002) Callus induction from Paulownia plant leaves and their plantlet regenerations. Sci Silv Sin 38: 29–35.

37. Hafner M, Landgraf P, Ludwig J, Rice A, Ojo T, et al. (2008) Identification of microRNAs and other small regulatory RNAs using cDNA library sequencing. Methods 44: 3–12.

38. Li R, Li Y, Kristiansen K, Wang J (2008) SOAP: short oligonucleotide alignment program. Bioinformatics 24: 713–714.

39. Zuker M (2003) Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res 31: 3406–3415.

40. Meyers BC, Axtell MJ, Bartel B, Bartel DP, Baulcombe D, et al. (2008) Criteria for annotation of plant MicroRNAs. Plant Cell 20: 3186–3190.

41. Audic S, Claverie JM (1997) The significance of digital gene expression profiles. Genome Res 7: 986–995.

42. Man MZ, Wang X, Wang Y (2000) POWER_SAGE: comparing statistical tests for SAGE experiments. Bioinformatics 16: 953–959.

43. Chen C, Ridzon DA, Broomer AJ, Zhou Z, Lee DH, et al. (2005) Real-time quantification of microRNAs by stem-loop RT-PCR. Nucleic Acids Res 33: e179.

44. Livak KJ, Schmittgen TD (2001) Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)). Methods 25: 402–408.

45. Wang K, Li M, Gao F, Li S, Zhu Y, et al. (2012) Identification of Conserved and Novel microRNAs from Liriodendron chinense Floral Tissues. PloS one 7: e44696.

46. Chen X (2009) Small RNAs and their roles in plant development. Ann Rev Cell Dev 25: 21–44.

47. Morin RD, O'Connor MD, Griffith M, Kuchenbauer F, Delaney A (2008) Application of massively parallel sequencing to microRNA profiling and discovery in human embryonic stem cells. Genome Res 18: 610–621.

48. Cifuentes M, Grandont L, Moore G, Chevre AM, Jenczewski E (2010) Genetic regulation of meiosis in polyploid species: new insights into an old question. New Phytol 186: 29–36.

49. Chen ZJ (2007) Genetic and epigenetic mechanisms for gene expression and phenotypic variation in plant polyploids. Annu Rev Plant Biol 58: 377–406.

50. Leitch AR, Leitch IJ (2008) Genomic plasticity and the diversity of polyploid plants. Science 320: 481–483.

51. Adams K (2007) Evolution of duplicate gene expression in polyploid and hybrid plants. J hered 98: 136–141.

52. Vashisht D, Nodine MD (2014) MicroRNA functions in plant embryos. Biochem Soc Trans 42: 352–357.

53. Ren Y, Chen L, Zhang Y, Kang X, Zhang Z, et al. (2012) Identification of novel and conserved Populus tomentosa microRNA as components of a response to water stress. Funct Integr Genomics 12: 327–339.

54. Zhao M, Ding H, Zhu J, Zhang F, Li W (2011) Involvement of miR169 in the nitrogen-starvation responses in Arabidopsis. New Phytol 190: 906–915.

55. Li B, Yin W, Xia X (2009) Identification of microRNAs and their targets from Populus euphratica. Biochem Bioph Res Co 388: 272–277.

56. Liu HH, Tian X, Li YJ, Wu CA, Zheng CC (2008) Microarray-based analysis of stress-regulated microRNAs in Arabidopsis thaliana. RNA 14: 836–843.

57. Mao W, Li Z, Xia X, Li Y, Yu J (2012) A combined approach of high-throughput sequencing and degradome analysis reveals tissue specific expression of microRNAs and their targets in cucumber. PloS one 7: e33040.

58. Chi X, Yang Q, Chen X, Wang J, Pan L, et al. (2011) Identification and characterization of microRNAs from peanut (Arachis hypogaea L.) by high-throughput sequencing. PLoS One 6: e27530.

59. Zhang X, Fan G, Zhao Z, Cao X, Zhao G, et al. (2013) Analysis of Diploid and Its Autotetraploid Paulownia tomentosa× P. fortunei with AFLP and MSAP. Sci Silv Sin 10: 026.

60. Shen H, He H, Li J, Chen W, Wang X, et al. (2012) Genome-wide analysis of DNA methylation and gene expression changes in two Arabidopsis ecotypes and their reciprocal hybrids. Plant Cell 24: 875–892.

61. Salmon A, Ainouche ML, Wendel JF (2005) Genetic and epigenetic consequences of recent hybridization and polyploidy in Spartina (Poaceae). Mol Ecol 14: 1163–1175.

62. Ha M, Lu J, Tian L, Ramachandran V, Kasschau K D, et al. (2009) Small RNAs serve as a genetic buffer against genomic shock in Arabidopsis interspecific hybrids and allopolyploids. Proc Natl Acad Sci USA 106: 17835–17840.

63. An F, Chan M (2012) Transcriptome-wide characterization of miRNA-directed and non-miRNA-directed endonucleolytic cleavage using Degradome analysis under low ambient temperature in Phalaenopsis aphrodite subsp. formosana. Plant Cell Physiol 53: 1737–1750

64. Wei M, Wei H, Wu M, Song M, Zhang J, et al. (2013) Comparative expression profiling of miRNA during anther development in genetic male sterile and wild type cotton. BMC Plant Biol 13: 66.

65. Rhoades MW, Reinhart BJ, Lim LP, Burge CB, Bartel B, et al. (2002) Prediction of plant microRNA targets. Cell 110: 513–520.

66. Nagpal P, Ellis CM, Weber H, Ploense SE, Barkawi LS, et al. (2005) Auxin response factors ARF6 and ARF8 promote jasmonic acid production and flower maturation. Development 132: 4107–4118.

67. Kwak P, Wang Q, Chen X, Qiu C, Yang Z (2009) Enrichment of a set of microRNAs during the cotton fiber development. BMC Genomics 10: 457.

68. Weir I, Lu J, Cook H, Causier B, Schwarz-Sommer Z, et al. (2004) CUPULIFORMIS establishes lateral organ boundaries in Antirrhinum. Development 131: 915–922.

69. Kosugi S, Ohashi Y (2002) DNA binding and dimerization specificity and potential targets for the TCP protein family. Plant J 30: 337–348.

70. Cubas P (2004) Floral zygomorphy, the recurring evolution of a successful trait. Bioessays 26: 1175–1184.

71. Koyama T, Furutani M, Tasaka M, Ohme-Takagi M (2007) TCP transcription factors control the morphology of shoot lateral organs via negative regulation of the expression of boundary-specific genes in Arabidopsis. Plant Cell 19: 473–484.

72. Sabatini S, Beis D, Wolkenfelt H, Murfett J, Guilfoyle T, et al. (1999) An auxin-dependent distal organizer of pattern and polarity in the Arabidopsis root. Cell 99: 463–472.

73. Sánchez C, Vielba JM, Ferro E, Covelo G, Solé A, et al. (2007) Two SCARECROW-LIKE genes are induced in response to exogenous auxin in rooting-competent cuttings of distantly related forest species. Tree Physiol 27: 1459–1470.

74. Wysocka-Diller JW, Helariutta Y, Fukaki H, Malamy JE, Benfey PN (2000) Molecular analysis of SCARECROW function reveals a radial patterning mechanism common to root and shoot. Development 127: 595–603.

75. Bouwmeester K, Govers F (2009) Arabidopsis L-type lectin receptor kinases: phylogeny, classification, and expression profiles. J Exp Bot 60: 4383–4396.

76. Singh P, Zimmerli L (2013) Lectin receptor kinases in plant innate immunity. Front Plant Sci 4: 124.

77. Li J, Jia D, Chen X (2001) HUA1, a regulator of stamen and carpel identities in Arabidopsis, codes for a nuclear RNA binding protein. Plant Cell 13: 2269–2281.

78. Berg JM, Shi Y (1996) The galvanization of biology: a growing appreciation for the roles of zinc. Science 271: 1081–1085.

79. Wang D, Guo Y, Wu C, Yang G, Li Y, et al. (2008) Genome-wide analysis of CCCH zinc finger family in Arabidopsis and rice. BMC Genomics 9: 44.

80. Peng X, Zhao Y, Cao J, Zhang W, Jiang H, et al. (2012) CCCH-type zinc finger family in maize: genome-wide identification, classification and expression profiling under abscisic acid and drought treatments. PLoS One 7: e40120.