Skip to main content

Identification of soybean seed developmental stage-specific and tissue-specific miRNA targets by degradome sequencing



MicroRNAs (miRNAs) regulate the expression of target genes by mediating gene silencing in both plants and animals. The miRNA targets have been extensively investigated in Arabidopsis and rice using computational prediction, experimental validation by overexpression in transgenic plants, and by degradome or PARE (parallel analysis of RNA ends) sequencing. However, miRNA targets mostly remain unknown in soybean (Glycine max). More specifically miRNA mediated gene regulation at different seed developmental stages in soybean is largely unexplored. In order to dissect miRNA guided gene regulation in soybean developing seeds, we performed a transcriptome-wide experimental method using degradome sequencing to directly detect cleaved miRNA targets.


In this study, degradome libraries were separately prepared from immature soybean cotyledons representing three stages of development and from seed coats of two stages. Sequencing and analysis of 10 to 40 million reads from each library resulted in identification of 183 different targets for 53 known soybean miRNAs. Among these, some were found only in the cotyledons representing cleavage by 25 miRNAs and others were found only in the seed coats reflecting cleavage by 12 miRNAs. A large number of targets for 16 miRNAs families were identified in both tissues irrespective of the stage. Interestingly, we identified more miRNA targets in the desiccating cotyledons of late seed maturation than in immature seed. We validated four different auxin response factor genes as targets for gma-miR160 via RNA ligase mediated 5’ rapid amplification of cDNA ends (RLM-5’RACE). Gene Ontology (GO) analysis indicated the involvement of miRNA target genes in various cellular processes during seed development.


The miRNA targets in both the cotyledons and seed coats of several stages of soybean seed development have been elucidated by experimental evidence from comprehensive, high throughput sequencing of the enriched fragments resulting from miRNA-guided cleavage of messenger RNAs. Nearly 50% of the miRNA targets were transcription factors in pathways that are likely important in setting or maintaining the developmental program leading to high quality soybean seeds that are one of the dominant sources of protein and oil in world markets.


MicroRNAs (miRNAs) are endogenous noncoding small RNAs which play significant roles in the regulation of gene expression. Post-transcriptional gene regulation by miRNAs constitutes one of the most conserved and well characterized gene regulatory mechanisms. It is important for growth, development, stress responses and numerous other biological processes in eukaryotes [15]. Therefore, identification of miRNAs and their targets in diverse species has been a major focus in recent years [6, 7]. In higher plants, miRNAs play significant roles in different developmental stages by regulating gene expression at transcriptional and post-transcriptional levels [812].

Most plant miRNAs facilitate the degradation of their mRNA targets by slicing precisely between the tenth and eleventh nucleotides (nt) from the 5’ end of the miRNA. As a result, the 3’ fragment of the target mRNA possesses a monophosphate at its 5’ end. This important property has been used to validate miRNA targets [8]. Isolation of such fragments is one of the critical steps for validating miRNA guided cleavage of target mRNA. A major limitation of this procedure is that every single predicted gene has to be verified separately. So, one-at-a-time isolation of cleaved RNA fragments is laborious, time-consuming and expensive. Recently, high-throughput sequencing methods, known as degradome analysis or PARE (parallel analysis of RNA ends) that can globally identify small RNA targets have been developed to overcome such limitations [13, 14].

Soybean is one of the most important crops cultivated all over the world. It is a good source of vegetable protein and oil. However, the role of miRNAs in soybean seed development is mostly unknown. So it is important to identify the seed developmental stage-specific and tissues-specific miRNAs and their potential target genes. Identification of the consequences of miRNA-guided target degradation that occurs in a developmental and tissue specific manner could help to elucidate how lipid and protein metabolic pathways operated during seed development. The soybean genome (cv. Williams82) was decoded a year ago [15], and this information has accelerated molecular research on soybeans. Although many soybean miRNAs were identified in previous research [1620], the number of miRNAs known in soybean is still very small and considerably lower than that in Arabidopsis or rice. High-throughput sequencing technologies such as massively parallel signature sequencing (MPSS), 454 and sequencing-by-synthesis (SBS) have enabled the identification of miRNAs in soybean. The extent of miRNA-directed post-transcriptional gene regulation in any organism can only be fully realized by identifying not only the miRNA component but also the set of their RNA targets.

Recently, miRNA targets have been reported for one of the many stages of soybean seed development, namely very early at 15 days after flowering, and without dissection of the maternal seed coats from the cotyledons which develop from the zygote [21]. To comprehensively investigate small RNA targets and provide basic information for further understanding of the miRNA-mediated post-transcriptional regulation during different soybean seed developmental stages, we constructed five separate degradome libraries derived from seed coats and cotyledons of different developmental stages representing the early, mid, and late maturation stages of seed development. The libraries were sequenced using SBS sequencing technology. The degradome dataset for the five different libraries was computationally analyzed. The majority of these reads mapped to the soybean transcriptome. A total of 183 target genes were confirmed as miRNA targets, which included both conserved and non-conserved miRNAs. Additionally, we have identified targets for 25 cotyledon-specific miRNAs, as well as 12 miRNAs and their potential targets found only in the seed coats. We found 16 miRNA families and their large number of targets that are found in both tissues. Moreover, we have validated Auxin Response Factors (ARFs) to be targets of gma-miR160, as verified by RNA ligase-mediated 5’ rapid amplification of cDNA ends (RLM 5’-RACE). The identification of developmental stage-specific and tissue-specific miRNA targets including many transcription factors advance our understanding of the miRNA-mediated post-transcriptional gene regulation during soybean seed development.


Seed developmental stage-specific library construction, sequencing and sequence analysis

In higher plants, most miRNAs regulate their targets via cleavage, which normally occurs between the tenth and eleventh nucleotides of the complementary region between the miRNA and the mRNA target [22]. The 3’ cleavage fragments contain both a free 5’ monophosphate and a 3’ polyA tail. So, these cleavage products can be successfully ligated with RNA ligase, whereas full length cDNAs with a 5’ cap structure or other RNAs lacking the 5’ monophosphate group are not compatible for ligation [8] and thus will be unavailable for subsequent amplification and sequencing reactions.

Five different degradome libraries, which capture the cleaved mRNAs, were constructed from cotyledons and seed coats from different seed developmental stages. These represented early maturation seed (25–50 mg fresh weight, green seed) and mid-maturation (100–200 mg fresh weight, green seed), the stages when the biosynthetic capacity of the seed is maximal and proteins and oils are accumulated at a high rate. In addition, we constructed a library from the yellow cotyledons (300–400 mg fresh weight) that are undergoing dehydration and desiccation.

SBS sequencing of these libraries produced raw reads from 10 million to 45 million (Table 1). After removal of low quality sequences and adapter removal, 95% of these reads lengths were 20 or 21 nt in length as expected from the cloning procedure. More than 97% of reads mapped to the soybean genome available at the Phytozome data base [23]. We also used the computationally predicted cDNA transcripts from the soybean genome sequence consisting of 78,773 high and low confidence gene models (Glyma models) for mapping degradome reads and found that more than 95% of reads matched to the Glyma models. These data indicate our degradome libraries to be of high quality and efficiency in recovering degraded mRNA targets that should contain the sequence profile resulting from miRNA directed cleavage.

Table 1 Analysis of degradome reads from five different libraries matched to the soybean genome

Systematic identification of miRNA targets in soybean

Systematic identification of miRNA targets was accomplished using previously described methods by analyzing the 20 and 21 nt reads with the CleaveLand pipeline for miRNA target identification [24] using all Glycine max miRNAs from miRBase [25, 26]. The identified targets were grouped into five categories by the program based on the relative abundance of the number of reads mapping to the predicted miRNA target site relative to other sites in the gene model (see Methods for details). Those in category 0 clearly have the majority of tags located at the miRNA-guided cleavage site.

The identified miRNA targets using degradome sequencing are presented in the form of target plots (t-plots) that plot the abundance of the signatures relative to their position in the transcript [14]. Representative t-plots are shown, one from each of four different degradome libraries such as cotyledon (25–50 mg), seed coat (25–50 mg), cotyledon (100–200 mg) and seed coat (100–200 mg) (Figure 1). In each of the four Glyma models, a clear peak for the absolute number of tags is found at the predicted cleavage site for gma-miR159, gma-miR4406, gma-miR167, or gma-miR164.

Figure 1
figure 1

Target plots (t-plots) of identified miRNA targets using degradome sequencing. Representative t-plots (a-d) are shown, one from each of four different libraries of the cotyledon (25–50 mg fresh weight range), seed coat (25–50 mg), cotyledon (100–200 mg) and seed coat (100–200 mg). Red arrows indicate signatures consistent with miRNA-directed cleavage. mRNA:miRNA alignments along with the detected cleavage frequencies (absolute numbers) are shown above the black arrow. The red colored italicized nucleotide on the target transcript from the 3’ end indicates the cleavage site detected in the degradome library.

The location of the cleavage site within the target gene is another important aspect of miRNA-mediated gene silencing. In soybean, the cleavage site of the miRNA was usually located in the CDS (coding sequence) of the target genes (Tables 23 and 4 and Additional file 1). Since the soybean genome at Phytozome [23] used computational predictions of gene models (known as Glyma models), some are likely deficient at the 5’ and 3’ UTRs (untranslated regions). Due to the some gene models being incomplete in the UTRs, there are likely other genes targeted by miRNA-guided cleavage in the UTR regions that may not be detected in our alignment analyses. In addition, miRNAs that function through translational repression, as opposed to cleavage of the target mRNA, will also not be identified by degradome or PARE sequencing techniques.

Table 2 Identified miRNA targets found in both soybean seed coat and cotyledon
Table 3 Potential miRNA targets found only in cotyledons at different seed developmental stages
Table 4 Potential miRNAs and their targets found only in seed coats at different developmental stages

The full complement of targets found in each of the five degradome libraries is presented in Additional file 1. In total, 183 targets representing 53 different miRNAs families were identified. Of those 133 targets were found representing the putative action of 16 different miRNAs in common between both tissues. Table 2 presents a subset of those that are found in at least one stage of development for both seed coats and cotyledons. The CleaveLand program predicts any gene family members that have a splice site matching the degradome data. Some miRNA family members residing at different genomic locations have very similar, if not identical mature miRNA sequences. Thus, the predictions from analysis of degradome data do not necessarily mean that the particular miRNA family member revealed from degradome data is the one expressed in that tissue. Direct sequencing of the small RNA population is required to verify the presence of a particular gene family member. Inspection of small RNA sequencing data from seed coats and cotyledons of Williams [27] shows the presence of various miRNA family members for gma-miR156, 159, 160, 164, 166, and 167, thus confirming that these miRNAs are present during seed development.

Identification of miRNA targets specific to either seed coat or cotyledons during seed development

Tissue specific miRNA and target identification is very important for understanding the regulation of gene expression in a spatial manner. In this study, we constructed cotyledon and seed coat libraries separately to identify miRNA targets both at younger (25–50 mg fresh weight seed) and older stages (100–200 mg and 300–400 mg fresh weight seed) of soybean seed development. Tissue specific siRNAs generated from a cluster of inverted repeat chalcone synthase (CHS) genes that downregulate CHS mRNAs and lead to lack of pigment on soybean seed coats have been described [27, 28], but very little is known about the miRNAs and their targets in developing seed tissues. We analyzed the degradome data from seed coats versus cotyledons and identified 25 miRNAs and their 32 different targets that were found only in the cotyledons and not the seed coats (Table 3). Likewise, 12 miRNAs and their 18 targets are associated with the seed coats only (Table 4).

Validation of miRNA targets

We report here that many targets were captured by the degradome analysis, which provided experimental evidence to support previous computational predictions. Because of its polyploid genome, many soybean genes are present in multiple copies. As a result, some of the reads align to multiple members of the same gene family. To further confirm the degradome data for some of the family members, a RLM-5’ RACE experiment was performed to examine which family members were targeted by the miRNA for degradation. For gma-miR160 in the cotyledon (100–200 mg) degradome library (Table 2 and Additional file 1), we have identified five targets annotated as Auxin Response Factors (ARFs). Four of the five, namely Glyma12g08110.1, Glyma12g29720.1, Glyma14g33730.1 and Glyma11g20490.1, were also verified by RLM-5’RACE to be subjected to cleavage guided by gma-miR160 (Figure 2).

Figure 2
figure 2

Validation of Auxin Response Factors (ARFs) regulated by gma-miR160 in cotyledon. Confirmed targets (a-d) for gma-miR160 are presented in the form of target plots (t-plot) and alignments. Absolute numbers of signature sequences are indicated in the t-plot. Red arrows indicate signatures consistent with miRNA-directed cleavage. The black arrows indicate a site verified by RLM 5’-RACE and detected cleavage frequencies (absolute numbers) are shown above the arrow. The cleavage site is shown as a red letter. Cleavage frequency as determined by gene-specific 5’-RACE at the indicated position is shown in parenthesis.

GO analysis of miRNA target genes in soybean seed developmental stages

The identified targets for miRNAs in the three cotyledon degradome libraries were classified by their gene ontology (GO) using the AgriGO toolkit [29] (Figure 3). Higher percentages of these targets were found to be involved in developmental, reproductive, and regulatory and metabolic processes with respect to their proportions within the GO classification of all soybean cDNAs. The same general pattern is found for the targets predicted with the seed coats ( Additional file 2). The enrichment of the genes involved in developmental and regulatory processes may be consistent with the fact that the degradome libraries were constructed from different stages of developing soybean seeds. For the developing seeds, it is of utmost important to accumulate proteins and lipids that are subsequently used as the source of energy and amino acids for the germinating seedling. The corresponding miRNAs may regulate the expression of these target genes during different seed developmental stages in soybean through affecting various transcription factors that induce or shut off specific metabolic networks during the course of seed development. Interestingly, we identified more miRNA targets in the cotyledons of late seed maturation than earlier stages with a total of 92 different targets in the 300–400 desiccating, yellow seeds compared to 60 and 53 total in the early and mid-maturation, immature green seed respectively.

Figure 3
figure 3

GO analysis of miRNA target genes identified in cotyledons in different soybean seed developmental stages. Green bars indicate the enrichment of miRNA targets in GO terms. Blue bars indicate the percentage of total annotated soybean genes mapping to GO terms. Only the predicted target genes for miRNAs identified by degradome sequencing were considered.


Regulation of gene expression by miRNAs has been comprehensively investigated in animals and plants [3, 4, 22, 30, 31]. In the case of higher plants, Arabidopsis and rice miRNA targets have been widely studied by high throughput sequencing [13, 3234]. Soybean is a polyploid crop plant having a complex and large genome compared to Arabidopsis and rice. The number of identified miRNAs and their potential targets in soybean is limited. To date, degradome sequencing has been reported for only one soybean tissue, namely the very young whole seed extracted 15 days after flowering from the cultivar Heinong44 [21]. In order to study the regulation of gene expression during soybean seed development, we constructed and sequenced five distinct degradome libraries using cotyledons and seed coats as a tissue source from different seed developmental stages of the cultivar Williams. The seed coat is primarily maternal tissue while the cotyledons represent the embryo of the next generation. Both tissues have the same genotype in inbred soybean lines.

After computational analysis using the Cleaveland pipeline [24], we identified a total of 183 potential targets of 53 miRNA families in five different seed developmental stage specific degradome libraries. Subsequent analysis and identification of cotyledon and seed coat specific miRNAs and their targets give us a better understanding of the regulation of gene expression in a spatial manner during soybean seed development at later stages of seed development in a widely used Maturity Group III cultivar “Williams”. The soybean genome sequence and predicted gene models used in this study are derived from a closely related isoline Williams82 [15]. Validation of four Auxin Response Factor (ARF) genes as targets of gma-miR160 indicates degradome sequencing as an efficient strategy to identify miRNA targets in plants.

Comparison of the data in our report with the first soybean degradome data reported by Song et al. [21] from a very early stage of soybean seed development in the cultivar Heinong44 revealed a number of miRNA targets in common among the data sets including many transcription factors such as ARF, MYB, TCP, NF-Y, Growth Regulatory Factor, HD-ZIP, PPR, SBP and NAC family protein. In contrast, we found some other miRNA targets such as Permease family protein, LRR domain containing protein, transmembrane protein 14 C, Serine-Threonine kinase, BRE expressed protein and transcription factor TFIID only in our degradome libraries that represented later stages of seed development of the cultivar Williams. Overall there were 65 Glyma model targets in common between both of the data sets, 80 found only in the Heinong44 data set, and 118 found only in the Williams data sets.

The 183 identified targets from our degradome analyses belonged to 57 different annotation groupings ( Additional file 1). Comparing those to the Song et al. paper, 11 annotation groups are in common, 18 are unique to the Heinong44 very young seed, and 46 are unique to the five libraries of the Williams data sets (Table 5). Many of the identified soybean miRNA targets belong to diverse gene families of transcription factors such as ARFs, MYBs, TCPs, NACs, HD-ZIPs and NF-Y subunits (Tables 23 and 4 and Additional file 1). Many of these transcription factors are known to regulate diverse aspects of plant growth and development. Patterning and outgrowth of lateral organs in plants depend on the expression of HD-ZIP transcription factors that specify adaxial/upper cell fate [35, 36]. We identified a number of HD-ZIP transcription factors as targets for gma-miR166 in our degradome libraries. MYB family members in rice, which are targeted by miR159, appear to play an important role in response to the presence of abscisic acid (ABA) during plant embryonic development, suggesting their roles in seed development [37, 38]. In this study, we detected a number of MYB family transcription factors regulated by gma-miR159 (Additional file 1. The gma-miRNA156 family members target sites in numerous proteins containing the Squamosa Promoter Binding (SBP) domain. SBP and SBP-LIKE (SPLs) proteins play multiple roles in plant development [39, 40]. In Arabidopsis, rice and in some other plants, miR156 regulates leaf development by targeting Squamosa-Promoter Binding protein-like (SBP) transcription factors [31, 39]. The identification of SBP as a target of gma-miR160 may indicate the additional level of regulation for SBP during soybean seed development.

Table 5 Comparison of miRNA target annotations between the Williams and Heinong44 data sets

Our analyses of the early, mid and late maturation developmental stages of soybean show a number of targets similar to those found by Song et al., [21] in the very young seeds of the cultivar Heinong44 including the SPB transcription factors. One notable difference was the absence in our degradome data of miRNA172 targets which include members of the AP2 transcription factor family. From inspection of sequenced small RNA populations from the 50–75 mg seed coats and cotyledons of Williams [27], we find only a few occurrences of the miR172 family (less than 30 occurrences per million reads) while some family members of the miR156 family are highly abundant (99,000 per million) in the cotyledons. We speculate that miR172 and/or its targets may be more abundant in the very young seed used by the Song et al. group [21] and not prevalent in the mid-maturation seed that we have examined. In Arabidopsis, miR172 has been reported to be involved in the regulation of flowering time and floral development [40]. Alternatively, the AP2 factors may not be detected as targets in the degradome data if translational repression by miR172 is operative as has been shown in Arabidopsis flower development [41]. Nuclear Factor Y (NF-Y) was shown to control a variety of agronomically important traits, including drought tolerance, flowering time, and seed development [42]. We detected seven NF transcription factor YA subunit mRNAs specifically in seed coats that are targets of miR169 family members that occur in both seed coats and cotyledons ( Additional file 1). These targets may indicate some specific regulation of NF-YA transcription factors during soybean seed development.

To obtain a deeper understanding of soybean seed development, we investigated tissue specific miRNA target identification in the cotyledons and seed coats at different seed developmental stages. Based on the degradome data, we identified some miRNAs that may act differentially in the cotyledons versus the seed coats to degrade their targets (Tables 3 and 4). F-box proteins involved in auxin-stimulated protein degradation (TIR1-like) were among the identified targets specifically found in soybean cotyledon (Table 3). The LRR kinases have been reported to play important roles in plant development and brassinosteroid and ABA signaling [43]. We identified several LRR domain containing proteins as targets for gma-miR393, gma-miR1523 and gma-miR2109 (Tables 23 and 4). The presence of these miRNA targets implies their regulation during soybean seed development. As a storage organ, the soybean seed contains significant amounts of lipid and protein. Thus the regulation of energy metabolism is very important during seed development. We identified a number of targets in the soybean cotyledon such as NADP/FAD oxidoreductase, ribose-5-phosphate isomerase, GTPase activating proteins and ferredoxin related proteins which are related to energy metabolism. Both in the soybean cotyledon and seed coat, we found pentatricopeptide repeat (PPR) proteins as targets of miR1520 (Table 2 and Additional file 1) which regulates gene expression in the mitochondria and chloroplasts [44, 45]. Since we constructed our degradome libraries using cotyledons and seed coats from different seed developmental stages, we identified targets of miRNAs during a broad range of soybean seed development.

Auxin is an important phytohormone in higher plants. It acts as a key player in plant development [46]. As the transducer of auxin signaling, ARFs play vital roles in plant development, including shoot, root and flower formation [47, 48]. In Arabidopsis, miR160 and miR167 are involved in auxin signaling via regulation of ARF genes [1]. In rice, a number of ARF encoding genes have been identified which are regulated by osa-miR160 and osa-miR167, respectively [33, 34]. In our study, we identified a large number of ARF genes as targets for different miRNAs such as gma-miR160 and miR167. In the cotyledon (100–200 mg) degradome library ( Additional file 1), we identified five targets annotated as Auxin Response Factors (ARFs) for gma-miR160, and four of these, Glyma12g08110.1 Glyma12g29720.1 Glyma14g33730.1 and Glyma11g20490.1, were validated by RLM-5’RACE showing precise cleavage as expected (Figure 2). These results suggested that gma-miR160 could participate in auxin signaling via down-regulation of ARFs during soybean seed developmental stages. The cleaved mRNAs captured by the degradome procedure indicate that the levels of the ARF mRNA targets are likely to be decreased, but qRTPCR or RNA sequencing data would be needed to directly confirm the effect on mRNA levels for a particular ARF target gene.

In our degradome libraries, miRNA targets are involved in major transitions between each stage of seed development and transcription factors account for approximately half of these targets. Of 183 identified targets in our soybean degradome libraries, GO analysis for biological function indicates that these genes are mainly involved in developmental and metabolic processes (Figure 3 and Additional file 2). Enrichment of developmentally related genes as target miRNAs suggests the high level of regulation of gene expression during soybean seed development. The larger number of targets found in the 300–400 mg desiccating, yellow cotyledons of late maturation implies that post-transcriptional regulation by miRNAs may aid in shifting the developmental program of the immature soybean cotyledons from biosynthesis of storage reserves to a catabolic role in utilization of those reserves during seed germination and growth. The miRNA targets verified by degradome sequencing will provide useful information for understanding and revealing significant roles of miRNAs during soybean seed development.


Degradome sequencing is a valuable tool for the experimental confirmation of miRNA targets in higher plants. This method can reveal additional targets which are difficult to identify by computational prediction alone and confirm that the targets genes have been cleaved in specific tissues. Five degradome libraries from three different developmental stages identified 183 miRNA targets. Identification of soybean seed coat and cotyledon specific miRNA targets gives better understanding of tissue specific miRNA targets during seed development. In summary, the current study has confirmed a large set of targets that are subjected to miRNA guided degradation including many transcription factors and a surprisingly large number of targets in the late stages of cotyledon development. The data provides an avenue to explore more details about developmental stage specific miRNA targets that play critical roles in each of the important tissues during seed development.


Plant materials

Soybean (Glycine max cv. Williams) plants were grown in a greenhouse and seeds were collected at different developmental stages including early maturation, green 25–50 mg fresh weight seed, mid- maturation green 100–200 mg, and late maturation yellow 300–400 mg fresh weight seed. Immediately, cotyledons and seed coats were separated by dissecting whole seeds and then frozen in liquid nitrogen. Subsequently the tissue was freeze dried and stored at −80°C.

RNA extraction and degradome library construction

Total RNA was extracted from freeze dried cotyledons and seed coats using a modified McCarty method [49] using phenol-chloroform extraction and lithium chloride precipitation. Approximately 200 μg of total RNA was used to select poly (A) RNA using the Oligotex mRNA mini kit (Qiagen). The degradome libraries were constructed as previously described [14, 50]. Briefly, using T4 RNA ligase (Ambion), a 5’ RNA adapter (5’-GUUCAGAGUUCUACAGUCCGAC-3’, Dharmacon Inc.) was added to the cleavage products, which possess a free 5’-monophosphate at their 3’ termini. Then the ligated products were purified by Oligotex mRNA mini kit (Qiagen) and reverse transcribed using an oligo dT primer (5’-CGAGCACAGAATTAATACGACTTTTTTTTTTTTTTTTTTV-3’, Integrated DNA Technologies) via SuperScript II RT (Invitrogen). The generated cDNA was amplified for 6 cycles (94°C for 30 s, 60°C for 20 s, and 72°C for 3 min) with a pair of primers (forward, 5´-GTTCAGAGTTCTACAGTCCGAC-3´and reverse, 5´-CGAGCACAGAATTAATACGACT-3´, Integrated DNA Technologies) using Phusion Taq (New England Biolabs). The PCR products were digested with restriction enzyme Mme I (NEB). Next, a double stranded DNA adapter (top 5’-p- TGGAATTCTCGGGTGCCAAGG-3’, bottom 5’ –CCTT GGCACCCGAGAATTCCANN-3’, Integrated DNA Technologies) was ligated to the digested products using T4 DNA ligase [13]. Then the ligated products were selected based on size by running 10% polyacrylamide gel. The gel-purified products were used for the final PCR amplification (94°C for 30 s, 60°C for 20 s, and 72°C for 20 s) with primers (forward, 5’-AATGATACGGCGACCACCGAGATCTACACGTTCAGAGTTCTACAGTCCGA-3’, reverse, 5’-CAAGCAGAAG ACGGCATACGAGATCGTGATGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA-3’, Integrated DNA Technologies) for 20 cycles. Finally, PCR products were gel purified and subjected to SBS sequencing by the Illumina HiSeq2000 at the Keck Center, University of Illinois at Urbana-Champaign. The sequencing data of the five degradome libraries are available under NCBI-GEO [51] series accession no. GSE34433.

Initial processing and analysis of reads for different sequencing libraries

Degradome libraries were sequenced by the Illimuna HiSeq2000. The raw data were preprocessed to remove low quality reads and clip adapter sequences. Subsequently only 20–21 nt sequences with high quality scores were collected for analysis. The ultrafast Bowtie aligner [52] was used to map soybean degradome reads to the Phytozome Glycine max gene models [23]. The distinct reads that perfectly matched soybean transcript sequences remained. The CleaveLand pipeline [24] was used to find sliced miRNA targets using the Phytozome Glycine max gene models [23] and all Glycine max miRBase, release 17 [25, 26] containing 207 mature miRNA sequences as input. All alignments with scores up to 7 and no mismatches at the cleavage site (between the tenth and eleventh nucleotides) were considered candidate targets. This analysis was performed separately for all five libraries. The identified targets were grouped into five categories based on the relative abundance of the degradome signatures at the miRNA target sites as determined by the program that indicates the abundance of the fragments mapping at the predicted miRNA target site relative to the abundance of fragments found at other sites. In category 0, the most abundant tags are found at the predicted site of miRNA guided cleavage and there was only one maximum on the transcript. If there was more than one abundant tag, it is indicated as category 1. In category 2, the abundance of cleavage signatures was less than the maximum but higher than the median. It was grouped as category 3 when the abundance of cleavage signatures was equal to, or less than the median. Very low abundance signatures including only 1 read was categorized as category 4.


For mapping the cleavage site within the miRNA target, a modified procedure for RNA ligase mediated (RLM-5’ RACE) was performed using the FirstChoice RLM-RACE Kit (Ambion). Total RNA was isolated from cotyledons and seed coats separately at different soybean seed developmental stages. Poly (A) + RNAs were purified using Oligotex mRNA mini kit (Qiagen). A 5’ RNA adaptor (5'-GCUGAUGGCGAUG AAUGAACACUGCGUUUGCUGGCUUUGAUGAAA- 3', Ambion RACE kit) was ligated to approximately 100 ng of mRNA using T4 RNA ligase. The ligated mRNAs were then reverse transcribed using oligo(dT) primer via M-MLV reverse transcriptase (Ambion). Two rounds of 5’ RACE reactions were performed with two nested primers (outer, 5'-GCTGATGGCGATGAATGAACACTG-3’; inner, 5'-CGCGGATCCGAAC ACTGCGTTTGCTGGCTTTGATG-3’, Ambion) and two gene-specific outer and inner primers (Additional file 3). Subsequently PCR products were gel purified and sequenced.

Authors’ information

Department of Crop Sciences, University of Illinois, Urbana, Illinois, 61801, USA.


  1. Jones-Rhoades MW, Bartel DP, Bartel B: MicroRNAs and their regulatory roles in plants. Annu Rev Plant Biol. 2006, 57: 19-53. 10.1146/annurev.arplant.57.032905.105218.

    Article  CAS  PubMed  Google Scholar 

  2. Mallory AC, Vaucheret H: Functions of microRNAs and related small RNAs in plants. Nat Genet. 2006, 38 (Suppl): S31-S36.

    Article  CAS  PubMed  Google Scholar 

  3. Axtell MJ, Bowman JL: Evolution of plant microRNAs and their targets. Trends Plant Sci. 2008, 13: 343-349. 10.1016/j.tplants.2008.03.009.

    Article  CAS  PubMed  Google Scholar 

  4. Voinnet O: Origin, biogenesis, and activity of plant microRNAs. Cell. 2009, 136: 669-687. 10.1016/j.cell.2009.01.046.

    Article  CAS  PubMed  Google Scholar 

  5. Bartel DP: MicroRNAs: target recognition and regulatory functions. Cell. 2009, 136 (2): 215-233. 10.1016/j.cell.2009.01.002.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  6. Yang T, Xue L, An L: Functional diversity of miRNA in plants. Plant Science. 2007, 172: 423-432. 10.1016/j.plantsci.2006.10.009.

    Article  CAS  Google Scholar 

  7. Sunkar R, Jagadeeswaran G: In silico identification of conserved microRNAs in large number of diverse plant species. BMC Plant Biol. 2008, 8: 37-10.1186/1471-2229-8-37.

    PubMed Central  Article  PubMed  Google Scholar 

  8. Llave C, Xie Z, Kasschau KD, Carrington JC: Cleavage of Scarecrow-like mRNA targets directed by a class of Arabidopsis miRNA. Science. 2002, 297: 2053-2056. 10.1126/science.1076311.

    Article  CAS  PubMed  Google Scholar 

  9. Aukerman MJ, Sakai H: Regulation of flowering time and floral organ identity by a microRNA and its APETALA2-like target genes. Plant Cell. 2003, 15: 2730-2741. 10.1105/tpc.016238.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  10. Chen X: MicroRNA biogenesis and function in plants. FEBS Lett. 2005, 579 (26): 5923-5931. 10.1016/j.febslet.2005.07.071.

    Article  CAS  PubMed  Google Scholar 

  11. Chuck G, Candela H, Hake S: Big impacts by small RNAs in plant development. Curr Opin Plant Biol. 2009, 12 (1): 81-86. 10.1016/j.pbi.2008.09.008.

    Article  CAS  PubMed  Google Scholar 

  12. Lelandais-Briere C, Sorin C, Declerck M, Benslimane A, Crespi M, Hartmann C: Small RNA diversity in plants and its impact in development. Curr Genomics. 2010, 11 (1): 14-23. 10.2174/138920210790217918.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  13. Addo-Quaye C, Eshoo TW, Bartel DP, Axtell MJ: Endogenous siRNA and miRNA targets identified by sequencing of the Arabidopsis degradome. Curr Biol. 2008, 18: 758-762. 10.1016/j.cub.2008.04.042.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  14. German MA, Pillay M, Jeong DH, Hetawal A, Luo SJ, Janardhanan P, Kannan V, Rymarquis LA, Nobuta K, German R, Paoli ED, Lu C, Schroth G, Meyers BC, Green PJ: Global identification of microRNA-target RNA pairs by parallel analysis of RNA ends. Nat Biotechnol. 2008, 26: 941-946. 10.1038/nbt1417.

    Article  CAS  PubMed  Google Scholar 

  15. Schmutz J, Cannon SB, Schlueter J, Ma J, Mitros T, Nelson W, Hyten DL, et al: Genome sequence of the palaeopolyploid soybean. Nature. 2010, 463: 178-183. 10.1038/nature08670.

    Article  CAS  PubMed  Google Scholar 

  16. Subramanian S, Fu Y, Sunkar R, Barbazuk WB, Zhu JK, Yu O: New and nodulation-regulated microRNAs in soybean roots. BMC Genomics. 2008, 9: 160-10.1186/1471-2164-9-160.

    PubMed Central  Article  PubMed  Google Scholar 

  17. Wang YW, Li PC, Cao XF, Wang XJ, Zhang A, Li X: Identification and expression analysis of miRNAs from nitrogen-fixing soybean nodules. Biochem Biophys Res Commun. 2009, 378: 799-803. 10.1016/j.bbrc.2008.11.140.

    Article  CAS  PubMed  Google Scholar 

  18. Chen R, Hu Z, Zhang H: Identification of microRNAs in wild soybean (Glycine soja). J Integr Plant Biol. 2009, 51: 1071-1079. 10.1111/j.1744-7909.2009.00887.x.

    Article  CAS  PubMed  Google Scholar 

  19. Joshi T, Yan Z, Libault M, Jeong DH, Park S, Green PJ, Sherrier DJ, Farmer A, May G, Meyers BC, Xu D, Stacey G: Prediction of new miRNAs and associated target genes in Glycine max. BMC Bioinforma. 2010, 11: S14-

    Article  Google Scholar 

  20. Kulcheski FR, Oliveira LF, Molina LG, Almerão MP, Rodrigues FA, Marcolino J, Barbosa JF, et al: Identification of novel soybean microRNAs involved in abiotic and biotic stresses. BMC Genomics. 2011, 12: 307-10.1186/1471-2164-12-307.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  21. Song QX, Liu YF, Hu XY, Zhang WK, Ma B, Chen SY, Zhang JS: Identification of miRNAs and their target genes in developing soybean seeds by deep sequencing. BMC Plant Biology. 2011, 11: 5-10.1186/1471-2229-11-5.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  22. Bartel DP: MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004, 116: 281-297. 10.1016/S0092-8674(04)00045-5.

    Article  CAS  PubMed  Google Scholar 

  23. Joint Genome Institute/Phytozome/. , , []

  24. Addo-Quaye C, Miller W, Axtell MJ: CleaveLand: a pipeline for using degradome data to find cleaved small RNA targets. Bioinformatics. 2009, 25: 130-131. 10.1093/bioinformatics/btn604.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  25. Griffiths-Jones S, Saini HK, van Dongen S, Enright AJ: miRBase: tools for microRNA genomics. Nucleic Acids Res. 2008, 36: D154-D158. 10.1093/nar/gkn221.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  26. Griffiths-Jones S, Grocock RJ, van Dongen S, Bateman A, Enright AJ: miRBase: microRNA sequences, targets and gene nomenclature. Nucleic Acids Res. 2006, 34: D140-D144. 10.1093/nar/gkj112.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  27. Tuteja JH, Zabala G, Varala K, Hudson M, Vodkin LO: Endogenous, Tissue-Specific Short Interfering RNAs Silence the Chalcone Synthase Gene Family in Glycine max Seed Coats. Plant Cell. 2009, 21: 3063-3077. 10.1105/tpc.109.069856.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  28. Kurauchia T, Kasaia A, Tougoub M, Senda M: Endogenous RNA interference of chalcone synthase genes in soybean: Formation of double-stranded RNA of GmIRCHS transcripts and structure of the 5’ and 3’ends of short interfering RNAs. Plant Physiol. 2011, 168: 1264-1270. 10.1016/j.jplph.2011.01.003.

    Article  Google Scholar 

  29. Du Z, Zhou X, Ling Y, Zhang ZH, Su Z: AgriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010, 38: W64-W70. 10.1093/nar/gkq310.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  30. Reinhart BJ, Weinstein EG, Rhoades MW, Bartel B, Bartel DP: MicroRNAs in plants. Genes Dev. 2002, 16: 1616-1626. 10.1101/gad.1004402.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  31. Rhoades MW, Reinhart BJ, Lim LP, Burge CB, Bartel B, Bartel DP: Prediction of plant microRNA targets. Cell. 2002, 110: 513-520. 10.1016/S0092-8674(02)00863-2.

    Article  CAS  PubMed  Google Scholar 

  32. Sunkar R, Zhou XF, Zheng Y, Zhang WX, Zhu JK: Identification of new and candidate miRNAs in rice by high throughput sequencing. BMC Plant Biol. 2008, 8: 25-10.1186/1471-2229-8-25.

    PubMed Central  Article  PubMed  Google Scholar 

  33. Li YF, Zheng Y, Addo-Quaye C, Zhang L, Saini A, Jagadeeswaran G, Axtell MJ, Zhang W, Sunkar R: Transcriptome-wide identification of microRNA targets in rice. Plant Journal. 2010, 62: 742-759. 10.1111/j.1365-313X.2010.04187.x.

    Article  CAS  PubMed  Google Scholar 

  34. Zhou M, Gu LF, Li PC, Song XW, Wei LY, Chen ZY, Cao XF: Degradome sequencing reveals endogenous small RNA targets in rice (Oryza sativa L. ssp. indica). Front Biol. 2010, 5: 67-90. 10.1007/s11515-010-0007-8.

    Article  CAS  Google Scholar 

  35. Emery JF, Floyd SK, Alvarez J, Eshed Y, Hawker NP, Izhaki A, Baum SF, Bowman JL: Radial Patterning of Arabidopsis Shoots by Class III HD-ZIP and KANADI Genes. Curr Biol. 2003, 13: 1768-1774. 10.1016/j.cub.2003.09.035.

    Article  CAS  PubMed  Google Scholar 

  36. Côté CL, Boileau F, Roy V, Ouellet M, Levasseur C, Morency MJ, Cooke JE, Séguin A, MacKay JJ: Gene family structure, expression and functional analysis of HD-Zip III genes in angiosperm and gymnosperm forest trees. BMC Plant Biology. 2010, 10: 273-10.1186/1471-2229-10-273.

    PubMed Central  Article  PubMed  Google Scholar 

  37. Aldridge CD, Probert RJ: Seed development, the accumulation of abscisic acid and desiccation tolerance in the aquatic grasses Porteresia coarctata (Roxb.) Tateoka and Oryza sativa L. Seed Science Research. 1993, 3: 97-103.

    Article  CAS  Google Scholar 

  38. Reyes J, Chua NH: ABA induction of miR159 controls transcript levels of two MYB factors during Arabidopsis seed germination. Plant J. 2007, 49: 592-606. 10.1111/j.1365-313X.2006.02980.x.

    Article  CAS  PubMed  Google Scholar 

  39. Zhang B, Pan X, Stellwag EJ: Identification of soybean microRNAs and their targets. Planta. 2008, 229: 161-182. 10.1007/s00425-008-0818-x.

    Article  CAS  PubMed  Google Scholar 

  40. Martin RC, Asalina M, Liu PP, Kristof JR, Coppersmith JL, Pluskota WE, Bassel GW, et al: The microRNA156 and microRNA172 gene regulation cascades at post-germinative stages in Arabidopsis. Seed Science Research. 2010, 20: 79-87. 10.1017/S0960258510000085.

    Article  CAS  Google Scholar 

  41. Chen X: A microRNA as a translational repressor of APETALA2 in Arabidopsis flower development. Science. 2004, 303: 2022-2025. 10.1126/science.1088060.

    Article  CAS  PubMed  Google Scholar 

  42. Cao S, Kumimoto RW, Siriwardana CL, Risinger JR, Holt BF: Identification and Characterization of NF-Y Transcription Factor Families in the Monocot Model Plant Brachypodium distachyon. PLoS One. 2011, 6: e21805-10.1371/journal.pone.0021805.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  43. Zhang XS, Choi JH, Heinz J, Chetty CS: Domain-Specific Positive Selection Contributes to the Evolution of Arabidopsis Leucine-Rich Repeat Receptor-Like Kinase (LRR RLK) Genes. J Mol Evol. 2006, 63: 612-621. 10.1007/s00239-005-0187-z.

    Article  CAS  PubMed  Google Scholar 

  44. Andrés C, Lurin C, Small LD: The multifarious roles of PPR proteins in plant mitochondrial gene expression. Physiol Plant. 2007, 129: 14-22. 10.1111/j.1399-3054.2006.00766.x.

    Article  Google Scholar 

  45. Schmitz-Linneweber C, Small I: Pentatricopeptide repeat proteins: a socket set for organelle gene expression. Trends Plant Sci. 2008, 13: 663-670. 10.1016/j.tplants.2008.10.001.

    Article  CAS  PubMed  Google Scholar 

  46. Rubio-Somoza I, Cuperus JT, Weigel D, Carrington JC: Regulation and functional specialization of small RNA–target nodes during plant development. Current Opinion in Plant Biology. 2009, 12: 622-627. 10.1016/j.pbi.2009.07.003.

    Article  CAS  PubMed  Google Scholar 

  47. Wang JW, Wang LJ, Mao YB, Cai WJ, Xue HW, Chen XY: Control of root cap formation by MicroRNA-targeted auxin response factors in Arabidopsis. Plant Cell. 2005, 17: 2204-2216. 10.1105/tpc.105.033076.

    PubMed Central  Article  CAS  PubMed  Google Scholar 

  48. Wang D, Pei K, Fu Y, Sun Z, Li S, Liu H, Tang K, Han B, Tao Y: Genome-wide analysis of the auxin response factors (ARF) gene family in rice (Oryza sativa). Gene. 2007, 394: 13-24. 10.1016/j.gene.2007.01.006.

    Article  CAS  PubMed  Google Scholar 

  49. McCarty DR: A simple method for extractions of RNA from maize tissue. Maize Genetics Cooperation News Letter. 1986, 60: 61-

    Google Scholar 

  50. German MA, Luo SJ, Schroth G, Meyers BC, Green PJ: Construction of Parallel Analysis of RNA Ends (PARE) libraries for the study of cleaved miRNA targets and the RNA degradome. Nat Protoc. 2009, 4: 356-362. 10.1038/nprot.2009.8.

    Article  CAS  PubMed  Google Scholar 

  51. Gene Expression Omnibus, a genomics repository database. , , []

  52. Langmead B, Trapnell C, Pop M, Salzberg SL: Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009, 10: R25-10.1186/gb-2009-10-3-r25.

    PubMed Central  Article  PubMed  Google Scholar 

Download references


We are grateful to Sean Bloomfield and Cameron Lowe for help with data analysis. We would like to extend our thanks to Sarah Jones for reading manuscript and making valuable comments. The research was funded by support from the United Soybean Board, the Illinois Soybean Association, and the USDA.

Author information

Authors and Affiliations


Corresponding author

Correspondence to Lila Vodkin.

Additional information

Competing interests

The authors declare that they have no competing interests.

Authors’ contributions

MS designed experiments, performed RNA extractions, constructed degradome libraries, analyzed degradome data and drafted the manuscript. LOV designed initial approach, obtained funding, led and coordinated the project, and edited the manuscript. All authors have read and approved the final manuscript.

Electronic supplementary material


Additional file 1: Complete list of identified miRNA targets identified in soybean seed coats and cotyledons during seed development. (XLS 154 KB)


Additional file 2: GO analysis of miRNA targets identified in seed coats in different soybean seed developmental stages. (DOC 66 KB)

Additional file 3: Gene Specific primer sequences for RLM-5’RACE. (XLS 21 KB)

Authors’ original submitted files for images

Below are the links to the authors’ original submitted files for images.

Authors’ original file for figure 1

Authors’ original file for figure 2

Authors’ original file for figure 3

Rights and permissions

This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Reprints and Permissions

About this article

Cite this article

Shamimuzzaman, M., Vodkin, L. Identification of soybean seed developmental stage-specific and tissue-specific miRNA targets by degradome sequencing. BMC Genomics 13, 310 (2012).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI:


  • Seed Coat
  • miRNA Target
  • Soybean Seed
  • Auxin Response Factor
  • Massively Parallel Signature Sequencing