Comprehensive analysis of RNA-seq data reveals the complexity of the transcriptome in Brassica rapa
© Tong et al.; licensee BioMed Central Ltd. 2013
Received: 22 April 2013
Accepted: 4 October 2013
Published: 7 October 2013
The species Brassica rapa (2n=20, AA) is an important vegetable and oilseed crop, and serves as an excellent model for genomic and evolutionary research in Brassica species. With the availability of whole genome sequence of B. rapa, it is essential to further determine the activity of all functional elements of the B. rapa genome and explore the transcriptome on a genome-wide scale. Here, RNA-seq data was employed to provide a genome-wide transcriptional landscape and characterization of the annotated and novel transcripts and alternative splicing events across tissues.
RNA-seq reads were generated using the Illumina platform from six different tissues (root, stem, leaf, flower, silique and callus) of the B. rapa accession Chiifu-401-42, the same line used for whole genome sequencing. First, these data detected the widespread transcription of the B. rapa genome, leading to the identification of numerous novel transcripts and definition of 5'/3' UTRs of known genes. Second, 78.8% of the total annotated genes were detected as expressed and 45.8% were constitutively expressed across all tissues. We further defined several groups of genes: housekeeping genes, tissue-specific expressed genes and co-expressed genes across tissues, which will serve as a valuable repository for future crop functional genomics research. Third, alternative splicing (AS) is estimated to occur in more than 29.4% of intron-containing B. rapa genes, and 65% of them were commonly detected in more than two tissues. Interestingly, genes with high rate of AS were over-represented in GO categories relating to transcriptional regulation and signal transduction, suggesting potential importance of AS for playing regulatory role in these genes. Further, we observed that intron retention (IR) is predominant in the AS events and seems to preferentially occurred in genes with short introns.
The high-resolution RNA-seq analysis provides a global transcriptional landscape as a complement to the B. rapa genome sequence, which will advance our understanding of the dynamics and complexity of the B. rapa transcriptome. The atlas of gene expression in different tissues will be useful for accelerating research on functional genomics and genome evolution in Brassica species.
KeywordsBrassica rapa RNA-seq Alternative splicing Transcriptome
The species Brassica rapa (2n=20, AA) includes several subspecies providing human nutrition in the form of leafy, root and stem vegetables and edible oil. It also represents the origin of the Brassica 'A’ genome and contributes to other cultivated oilseed crops of Brassica allopolyploids: B. napus (AACC) and B. juncea (AABB) . Therefore, B. rapa has great potential as a model for genomic and evolutionary research in Brassica species. Over the past decade, a growing number of genomic resources for B. rapa have become available [2–5], in particular the whole genome sequence of B. rapa accession Chiifu-401-42, an inbred Chinese cabbage line, which was published in late 2011 .
With the availability of the B. rapa genome sequence , it is essential to further identify and determine the activity of all functional elements on a genome-wide scale. For this purpose, a comprehensive analysis of the transcriptome is required to reveal potentially transcribed regions and understand expression patterns of the entire B. rapa gene model sets across tissues. High-throughput RNA-seq technology represents a powerful and cost-efficient tool for transcription profiling . First, it can detect and quantify gene expression with digital measurements, and is especially sensitive for low-expressed genes . Second, it allows researchers to improve gene annotation by extending transcriptional boundaries of genes and permits the discovery of novel genes or transcripts . Third, the most attractive advantage is that RNA-seq can survey alternative splicing (AS) events at single nucleotide resolution [10–13]. In addition, RNA-seq data show a high level of reproducibility in both technical and biological replicates .
In this study, we generated RNA-seq data from six different tissues (root, stem, leaf, flower, silique and callus) of the B. rapa accession Chiifu-401-42, which was the reference line used for whole genome sequencing . By performing a comprehensive analysis of these RNA-seq data from various tissues, our goals were: (i) to identify actively transcribed regions in the B. rapa genome, including annotated coding sequences and novel transcribed regions (NTRs); (ii) to unveil the transcriptome dynamics of B. rapa gene models, and identify candidate genes that are conserved, specific, or show correlated expression across tissues; and (iii) study the AS of B. rapa genes and investigate whether 'AS-preferred’ genes are enriched in certain functional categories or are associated with particular gene structure, such as short introns.
Results and Discussion
Generation, mapping, and assessment of RNA-seq reads
Paired-end RNA-seq reads of 90 bp in length were generated from six major organs or tissues of B. rapa including root (two samples), stem, leaf (two samples), flower, silique and callus. The base quality of RNA-seq reads was checked and analyzed using FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) (in Additional file 1: Figure S1 and Figure S2). We used relatively stringent criteria for quality control by removing the pair-end reads containing Ns or those where the number of bases whose PHRED-like score was less than 20 exceeded 10% (see Methods, in Additional file 1: Figure S3, in Additional file 2: Table S1). For eliminating the errors from the mis-priming of primers in RNA-seq experiments, the first nine base pairs of reads were trimmed before read alignment. Using the software TopHat2  and allowing an at most 2-bp mismatch, we obtained over 219 million uniquely mapped reads for subsequent analyses (in Additional file 2: Table S2). Seventy-four percent of uniquely mapped reads were contiguously mapped back to the B. rapa v1.5 genome sequences, and the others were splice-mapped to, and spanned potential splice junctions of the genome (in Additional file 1: Figure S4, in Additional file 2: Table S2). The uniquely mapped reads covered over 80 million bases of the total genomic bases, with an average read-depth of 26–75 per base across different samples (in Additional file 2: Table S2). Around 81% of uniquely mapped reads confirmed the expression of ~40 Mb out of a total of 50 Mb of annotated B. rapa coding sequence (CDS). The remaining ~16% and ~3% of uniquely mapped reads were aligned to another ~40 million genomic bases in intergenic and intronic regions, respectively (in Additional file 2: Table S2), implying that large amounts of actively transcribed regions remain unidentified in the B. rapa genome. The average read depth of coverage in CDS, introns and intergenic regions are shown in Additional file 1: Figure S4. As expected, the read depth in intergenic regions and introns was lower than that in exons, which correlated with the subsequent observation in this study that novel transcripts discovered from these regions are commonly expressed at low levels. The read depth seemed to be distributed relatively evenly along the whole body of genes (in Additional file 1: Figure S5 and Figure S6), reflecting no introduction of obvious bias during randomly primed reverse transcription and subsequent RNA sequencing. The two samples of root and leaf tissues from different individual plants both showed a high degree of concordance in gene expression of log2-transformed fragments per kilobase of transcript per million fragments mapped (FPKM) values (R-squared of 0.94 for both comparisons, in Additional file 1: Figure S7).
RNA-seq data improved the annotation of the B. rapa genome
The transcribed regions/units were constructed independently in individual tissues and then merged into a final dataset containing 37,002 transcripts (in Additional file 2: Table S3). 33,598 transcripts from this dataset completely (43%) or partially (57%) overlapped with 30,322 annotated gene models in B. rapa, and the remaining transcripts results in the identification of 3,450 novel transcribed regions (NTRs) in the B. rapa genome (in Additional file 2: Table S4).
Gene ontology (GO) categories with significantly longer or shorter UTRs were obtained by comparing the lengths of all UTRs of each GO category with total UTR length (Wilcoxon Rank Sum test, P < 0.05). The transcripts/genes related to regulatory signals and binding domains seem to have long 5′ and 3′-UTRs, as exemplified in the transcripts/genes encoding protein kinases, membrane proteins, and signal transducers (Figure 1B). By contrast, ribosomal genes seem to have short 5′-UTRs. This observation is consistent with previous reports in other organisms [11, 17]. In addition, we scanned a small upstream ORF (uORF) of start codons, and found that 4,147 (~29% of the total) genes contain a uORF (10–50 nt) in the 5′-UTR region (in Additional file 2: Table S6). The importance and mechanism of uORFs as regulatory factors for gene expression have been revealed clearly in recent years [18, 19]; therefore, further analyses of these novel UTRs should greatly aid research towards identifying all known regulatory elements in the B. rapa genome.
Discovery of novel transcripts
We identified 3,450 NTRs that are not linked to any annotated gene models in B. rapa (in Additional file 2: Table S4). Over 41.2% (1,422) of the NTRs have an ORF greater than 100 amino acids as predicted by Augustus , and 802 of them contain pfam protein domains . Additionally, 166 NTRs have sequence similarity (greater than 100-nucleotide match and 60% identity) to transposable elements (TEs) in Repbase (v16.10), reflecting substantial TE activity in the B. rapa genome. Only ~20% of them are supported by expressed sequence tags (ESTs) (using blat, more than 95% match in identity and coverage), despite using 214,106 ESTs from the National Center for Biotechnology Information (NCBI) database. To validate the predicted novel transcripts, we randomly selected six of them for reverse transcription polymerase chain reaction (RT-PCR) analysis. All six were confirmed to have transcriptional activity in at least one tissue (in Additional file 1: Figure S8).
These NTRs were difficult to discover previously because of their low-level and organ-specific expression; however, high-throughput sequencing technology and individual organ libraries permitted their discovery. Further analysis showed that highly expressed and poorly expressed novel transcripts were relatively equally detected in silique and callus tissues, unlike in other tissues, where poorly expressed novel transcripts were preferentially discovered (in Additional file 1: Figure S9). We speculated that because the current B. rapa genomic annotation is mainly based on ESTs collected from normal tissues, transcripts preferentially expressed in these two relatively specialized tissues could have remained unannotated. Thus, generating RNA-seq data from more diverse tissues from different developmental stages or environmental conditions will be helpful for improving the annotation of the B. rapa genomic sequence.
Tissue dynamics of B. rapa transcriptomes
For each gene model, it was considered expressed if lower boundary of the FPKM values for the 95% confidence interval of the abundance of gene was greater than zero (in Additional file 2: Table S7). 32,335 genes (78.8% of the total 41,020 annotated genes) were detected as expressed in at least one tissue, of which 18,876 genes were constitutively expressed in all tissues (Figure 3B). The clustering of log2-transformed FPKM values of constitutively expressed genes shows the transcriptomic relationship between tissues as a dendrogram: (((leaf, stem), root), (flower, silique)), callus) (Figure 3C).
We obtained a dataset containing 867 highly and stably expressed genes representing tissue-conserved expressed genes whose FPKM values are greater than 50 in every tissue (in Additional file 2: Table S8). This dataset contains potential useful housekeeping genes such as glyceraldehyde-3-phosphate dehydrogenase, actin, ubiquitin, tubulin and elongation factor 1-a, ribosomal protein, calmodulin, glutathione peroxidase and G-box regulating factor. For tissue-specific expressed genes, we observed some significantly meaningful gene groups in flower and silique tissues (in Additional file 2: Table S9). Many genes encoding lipid transfer proteins, oleosin, glycosyl hydrolase, lipase and transferases involved in protein amino acid glycosylation and phosphorylation, biosynthesis of flavonoids, biosynthesis of long chain fatty acids, and cellulose biosynthetic process were specifically expressed and functional in siliques. With respect to flower-specific gene expression, rapid alkalinization factor, pectate lyase, pollen Ole e 1 allergen and extensin family protein, pectinesterase family protein, arabinogalactan-protein involved pollen tube growth, invertase/pectin methylesterase inhibitor, pectate lyase, and S1 self-incompatibility protein-related were identified.
To identify highly correlated groups of genes (called 'modules’), we conducted weighted gene correlation network analysis using the R-package [25, 26]. After inspection of functional annotation of genes (inferred from the most similar gene in A. thaliana) in each module, we identified some meaningful co-expressed genes (Figure 3D, in Additional file 2: Table S10). First, four modules (ID: 1–4) abundant in ribosomal genes (40s or 60s subunit) were highly expressed across all tissues. Second, one module (ID: 5, 196 genes) is nearly exclusively committed to photosynthesis, with predominant expression of genes in leaf tissue encoding photosystem I and II subunit proteins, chlorophyll I and II binding proteins, rubisco subunit and activase, and other genes involved in chlorophyll biosynthesis, the photosynthetic electron transport chain, carbon utilization by fixation of carbon dioxide and light harvesting. We noticed that a similar group of co-expressed genes was reported in the maize leaf developmental transcriptome , implying that indispensable co-expression of these genes is required for leaf photosynthesis. Third, three modules (ID: 7–9) were significantly enriched for protein kinases and transcription factors (TFs) (146 genes), which were highly expressed in flowers.
Identification of splice junctions and alternative splicing (AS) events
A total of 156,516 unique splice junctions (read depth of each junction ≥ 3 in each sample) were detected using TopHat  (in Additional file 2: Table S11); 122,693 (78.4%) of them confirmed 74.1% of all known 165,564 exon-exon junctions in the B. rapa genome. Examination of the dinucleotides at the splicing borders showed that 98.4% of them had GT-AG splicing junctions (in Additional file 1: Figure S11).
Predicted alternative splicing (AS) events in B. rapa
Genes with AS event
AS events shared by more than two tissues
We further investigated the tissue distribution of predicted splice junctions (Figure 4B). We divided all 165,564 read-supported splice junctions into two classes: the 122,693 known or annotated genomic exon-exon splice junctions and the remaining 33,823 novel splice junctions. It could be assumed that the former represent the junctions used for normal transcription, whereas the latter (novel splice junctions) always reflect the existence of AS and produce novel splice variants of transcripts [30, 31]. We observed that each tissue dataset contributed between 50% and 78% of the read-supported known exon-exon splice junctions and the majority (90%) of them was detected to be shared by more than two tissues. This represents the normal and most frequent use of constitutive transcripts or isoforms in different tissues. For novel splice junctions, 65.6% of them were shared by more than two tissues and the remaining were detected in only one tissue representing tissue-specific or tissue-restricted splice variants (Figure 4B). It may be that AS events appear to be tissue-specific because they are rare, and more samples and replicates must been included to confirm the estimation for the tissue-specific splice junctions. However, a larger proportion of novel splice junctions are commonly detected in multiple tissues, indicating the prevalence of AS across tissues in the B. rapa genome.
Analysis of GO enrichment for AS genes, using all B. rapa genes as the background, revealed that GO categories enriched with AS genes were associated with signal transduction, regulation, response, binding and catalytic activity (Figure 4C, in Additional file 2: Table S13). It is worth noting that a GO category (GO: 0032318: regulation of Ras GTPase activity) and its paternal GO categories were significantly abundant in AS genes. Ras/small GTPases generally play a role in signal transduction as molecular switches for a variety of cellular signaling events , and AS may be essential for these genes to function in regulation and signal transduction.
Considering the possible AS effect on 'regulation’-related genes, we conducted a systematic survey of AS events that occurred in genes encoding TFs. First, based on annotated TF genes in A. thaliana (http://arabidopsis.med.ohio-state.edu/AtTFDB/), we identified 2,361 B. rapa genes and classified them into 57 different TF families through homology searching and domain identification (in Additional file 2: Table S14). Second, about 89% of TF genes showed expression in at least one tissue and 27.7% of those underwent AS. Subsequently, we examined the ratio of AS genes in each TF subfamily, and identified 11 TF subfamilies, such as ARF, AP2, MIKC, C3H and MYB-related, that have higher AS frequencies than others (Figure 4D). We found that previously reported AS events in TF genes, such as in CCA1/LHY-like MYB , R2R3-type MYB , MIKC-type MADS , bZIP , and NAC , almost all belong to these 'AS-preferred’ TF subfamilies. As suggested in metazoan , the above observation shows that alternative splicing may serve as an important and prevalent mechanism for these TF-encoding genes to form protein complexes that function in transcription in B. rapa.
Intriguingly, these 'AS-preferred’ genes were also those that were more likely to have been retained in mesopolyploid B. rapa during the diploidization process following a recent whole genome triplication that occurred in a Brassica ancestor ~13 million years ago . In other angiosperms TFs and signal transducer genes are preferentially retained as duplicates after whole genome duplication (WGD) events [45–48]. This was popularly explained by the gene balance hypothesis , which suggests that genes participate in the formation of macromolecular complexes, or are involved in more interaction or regulation are dosage sensitive and are thus expected to be retained more often after WGD . AS is able not only to remove functional domains to produce non-functional transcripts, thus regulating gene dosage [44, 51], but also causes subfunctionalization between WGD-derived paralogous genes [52, 53]. Thus, the observation that over-retained genes always have a high frequency of AS suggested that alternative splicing, as a fundamental aspect of the expression of many genes, may be a form of divergence of duplicated genes and thus contributed to their retention, although the possibility that AS and retention may are parallel and independent phenomena driven by a third selective factor, could not be excluded.
In conclusion, the RNA-seq data revealed the extent of transcription of the B. rapa genome, and simultaneously improved the existing gene annotation. The identification of numerous novel NTRs, commonly poorly or tissue-specifically expressed, led to the proposal that further transcriptome data from high-throughput platform will unveil previously unidentified functional regions of the B. rapa genome. A complete picture of the transcriptome across the major organs of B. rapa was provided, including candidate tissue-specific, tissue-conserved, or tissue-correlated expressed genes. The comprehensive analysis of AS and splice junctions increased our understanding of the complexity of the B. rapa transcriptome.
Organ collection and RNA sequencing
Six tissues of B. rapa accession Chiifu-401-42, including callus, root, stem, leaf, flower, and silique, were prepared for mRNA extractions. Plants were grown under greenhouse conditions at 22°C. Callus tissue was obtained from tissue culture. Root, stem and leaf tissues were collected from seven-week old plants. Two samples of root and leaf tissues were generated from different batches of plants. Flower tissue was obtained from blooming plants on the same day without the floral shoot. Silique tissue was generated from 15-day plant after pollination. We constructed an individual cDNA library with insert sizes of 200 bp for each sample, and sequenced them on the Illumina Hiseq 2000 platform. The libraries were sequenced for paired-end reads of 90 bp. All RNA-seq raw data from this study have been submitted to the NCBI Gene Expression Omnibus (http://www.ncbi.nlm.nih.gov/geo/) under accession no. GSE43245.
Mapping short reads and constructing transcripts
We used FastQC (http://www.bioinformatics.babraham.ac.uk/projects/fastqc/) to check and visualize the quality of the RNA-seq reads (in Additional file 1: Figure S1 and Figure S2). The average read quality was from 36 to 38 in different sample and the first nine base pairs of the reads, which showed unstable base composition given by the percentages of four different nucleotides (A, T, C, and G), were removed. We use the tools from the NGSQCToolkit (v2.3)  to remove the pair-end reads containing Ns or those where number of bases whose PHRED-like score was less than 20 exceeded 10% (in Additional file 2: Table S1). We calculated the sequence similarity of 10447 duplicated genes from the B. rapa triplicated genomes to evaluate the possible false reads for mapping (in Additional file 1: Figure S3). The B. rapa genome sequence (v1.5) was downloaded from the BRAD database (http://brassicadb.org/) for mapping and analysis We trimmed the first 9bp of filtered reads before read mapping using TopHat (v 2.0.9.) . The uniquely mapped reads were used for subsequent analysis, including transcripts construction, quantification of gene and transcript expression, and identification of splice junction and AS events. The transcripts were constructed and quantified the expression FPKM values of transcripts in each sample by Cufflinks (v2.1.0) and merged into a final comprehensive set of transcripts using Cuffmerge.
Tissue breadth of transcripts expression
In the formula, n is the number of surveyed B. rapa tissues. S (i, j) is the FPKM values of i gene in j tissue, and S (i, max) is the highest FPKM values of gene i in the n tissues.
Weighted Correlation Network Analysis (WGCNA)
Gene expression FPKM values were log2 transformed before being processed through the WGCNA R package . Co-expression analysis was performed to identify modules of highly correlated genes by the suggested protocol . We used the test to select the best parameters (softpower = 26). All other parameters were used with the default values.
Alternative splicing prediction
We used TopHat software  to identify splice junctions, then filtered junctions with read support of less than three. Each junction was searched for putative donor sites and acceptor sites inside intron regions. We predicted four types of AS events including exon skipping (ES), intron retention (IR), alternative 5′ splice site (A5SS) and alternative 3′ splice site (A3SS) according to the method described in Zhang et al. (2010) .
GO enrichment analysis
The Fisher's exact test (N ≤ 5) or Chi-square test (N > 5) was used to calculate the significance (P-value) for each GO category. The Benjamini-Hochberg false discovery rate (FDR) was performed to adjust the P-values . We used REVIGO (http://revigo.irb.hr/) to visualize enriched GO categories.
Identification of transcription factors in B. rapa
Putative transcription factors in B. rapa were identified by BLASTP search (E-value less than 1e-5), and redundant sequences were manually discarded. The TFs domains were scanned using HMMER 3.0 software (http://hmmer.janelia.org/) using the HMM profile, which was downloaded from Pfam and PlantTFDB databases (http://arabidopsis.med.ohio-state.edu/AtTFDB/).
Alternative 5' splice site donor
Alternative 3' splice site acceptor
- B. rapa:
Expressed sequence tags
Fragments per kilobase of exon per million fragments mapped
Novel transcribed regions
Open reading frame
Weighted gene correlation network analysis.
This research was supported by the National Natural Science Foundation of China (no. 31301039), the National High Technology Research and Development Program of China (863 Program, 2013AA102602), the National Basic Research and Development Program (973 Program, no. 2011CB109305), the Commonweal Specialized Research Fund of China Agriculture (201103016) and the Core Research Budget of the Non-profit Governmental Research Institution (no. 1610172011011).
- U N: Genome analysis in Brassica with special reference to the experimental formation of B. napus and peculiar mode of fertilication. Jap J Bot. 1935, 7: 389-452.
- Park JY, Koo DH, Hong CP, Lee SJ, Jeon JW, Lee SH, Yun PY, Park BS, Kim HR, Bang JW: Physical mapping and microsynteny of Brassica rapa ssp. pekinensis genome corresponding to a 222 kbp gene-rich region of Arabidopsis chromosome 4 and partially duplicated on chromosome 5. Mol Genet Genomics. 2005, 274 (6): 579-588. 10.1007/s00438-005-0041-4.View ArticlePubMed
- Mun JH, Kwon SJ, Yang TJ, Kim HS, Choi BS, Baek S, Kim JS, Jin M, Kim JA, Lim MH: The first generation of a BAC-based physical map of Brassica rapa. BMC Genomics. 2008, 9: 280-10.1186/1471-2164-9-280.PubMed CentralView ArticlePubMed
- Xiong Z, Kim JS, Pires JC: Integration of genetic, physical, and cytogenetic maps for Brassica rapa chromosome A7. Cytogenet Genome Res. 2010, 129 (1–3): 190-198.View ArticlePubMed
- Mun JH, Kwon SJ, Seol YJ, Kim JA, Jin M, Kim JS, Lim MH, Lee SI, Hong JK, Park TH: Sequence and structure of Brassica rapa chromosome A3. Genome Biol. 2010, 11 (9): R94-10.1186/gb-2010-11-9-r94.PubMed CentralView ArticlePubMed
- Wang X, Wang H, Wang J, Sun R, Wu J, Liu S, Bai Y, Mun JH, Bancroft I, Cheng F: The genome of the mesopolyploid crop species Brassica rapa. Nat Genet. 2011, 43 (10): 1035-1039. 10.1038/ng.919.View ArticlePubMed
- Wang Z, Gerstein M, Snyder M: RNA-Seq: a revolutionary tool for transcriptomics. Nat Rev Genet. 2009, 10 (1): 57-63. 10.1038/nrg2484.PubMed CentralView ArticlePubMed
- Mortazavi A, Williams BA, McCue K, Schaeffer L, Wold B: Mapping and quantifying mammalian transcriptomes by RNA-Seq. Nat Methods. 2008, 5 (7): 621-628. 10.1038/nmeth.1226.View ArticlePubMed
- Ozsolak F, Milos PM: RNA sequencing: advances, challenges and opportunities. Nat Rev Genet. 2011, 12 (2): 87-98. 10.1038/nrg2934.PubMed CentralView ArticlePubMed
- Sultan M, Schulz MH, Richard H, Magen A, Klingenhoff A, Scherf M, Seifert M, Borodina T, Soldatov A, Parkhomchuk D: A global view of gene activity and alternative splicing by deep sequencing of the human transcriptome. Science. 2008, 321 (5891): 956-960. 10.1126/science.1160342.View ArticlePubMed
- Wilhelm BT, Marguerat S, Watt S, Schubert F, Wood V, Goodhead I, Penkett CJ, Rogers J, Bahler J: Dynamic repertoire of a eukaryotic transcriptome surveyed at single-nucleotide resolution. Nature. 2008, 453 (7199): 1239-1243. 10.1038/nature07002.View ArticlePubMed
- Graveley BR, Brooks AN, Carlson JW, Duff MO, Landolin JM, Yang L, Artieri CG, van Baren MJ, Boley N, Booth BW: The developmental transcriptome of Drosophila melanogaster. Nature. 2011, 471 (7339): 473-479. 10.1038/nature09715.PubMed CentralView ArticlePubMed
- Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L: Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc. 2012, 7 (3): 562-578. 10.1038/nprot.2012.016.PubMed CentralView ArticlePubMed
- Marioni JC, Mason CE, Mane SM, Stephens M, Gilad Y: RNA-seq: an assessment of technical reproducibility and comparison with gene expression arrays. Genome Res. 2008, 18 (9): 1509-1517. 10.1101/gr.079558.108.PubMed CentralView ArticlePubMed
- Trapnell C, Pachter L, Salzberg SL: TopHat: discovering splice junctions with RNA-Seq. Bioinformatics. 2009, 25 (9): 1105-1111. 10.1093/bioinformatics/btp120.PubMed CentralView ArticlePubMed
- Alexandrov NN, Troukhan ME, Brover VV, Tatarinova T, Flavell RB, Feldmann KA: Features of Arabidopsis genes and genome discovered using full-length cDNAs. Plant Mol Biol. 2006, 60 (1): 69-85. 10.1007/s11103-005-2564-9.View ArticlePubMed
- Wang B, Guo G, Wang C, Lin Y, Wang X, Zhao M, Guo Y, He M, Zhang Y, Pan L: Survey of the transcriptome of Aspergillus oryzae via massively parallel mRNA sequencing. Nucleic Acids Res. 2010, 38 (15): 5075-5087. 10.1093/nar/gkq256.PubMed CentralView ArticlePubMed
- Jorgensen RA, Dorantes-Acosta AE: Conserved peptide upstream open reading frames are associated with regulatory genes in angiosperms. Front Plant Sci. 2012, 3: 191-PubMed CentralPubMed
- Somers J, Poyry T, Willis AE: A perspective on mammalian upstream open reading frame function. Int J Biochem Cell Biol. 2013, 45 (8): 1690-1700. 10.1016/j.biocel.2013.04.020.View ArticlePubMed
- Stanke M, Schoffmann O, Morgenstern B, Waack S: Gene prediction in eukaryotes with a generalized hidden Markov model that uses hints from external sources. BMC Bioinforma. 2006, 7: 62-10.1186/1471-2105-7-62.View Article
- Quevillon E, Silventoinen V, Pillai S, Harte N, Mulder N, Apweiler R, Lopez R: InterProScan: protein domains identifier. Nucleic Acids Res. 2005, 33: W116-W120. 10.1093/nar/gki442.PubMed CentralView ArticlePubMed
- Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L: Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010, 28 (5): 511-515. 10.1038/nbt.1621.PubMed CentralView ArticlePubMed
- Yanai I, Benjamin H, Shmoish M, Chalifa-Caspi V, Shklar M, Ophir R, Bar-Even A, Horn-Saban S, Safran M, Domany E: Genome-wide midrange transcription profiles reveal expression level relationships in human tissue specification. Bioinformatics. 2005, 21 (5): 650-659. 10.1093/bioinformatics/bti042.View ArticlePubMed
- Zhang G, Guo G, Hu X, Zhang Y, Li Q, Li R, Zhuang R, Lu Z, He Z, Fang X: Deep RNA sequencing at single base-pair resolution reveals high complexity of the rice transcriptome. Genome Res. 2010, 20 (5): 646-654. 10.1101/gr.100677.109.PubMed CentralView ArticlePubMed
- Zhang B, Horvath S: A general framework for weighted gene coexpression network analysis. Stat Appl Genet Mol Biol. 2005, 4: 17-
- Langfelder P, Horvath S: WGCNA: an R package for weighted correlation network analysis. BMC Bioinforma. 2008, 9: 559-10.1186/1471-2105-9-559.View Article
- Li P, Ponnala L, Gandotra N, Wang L, Si Y, Tausta SL, Kebrom TH, Provart N, Patel R, Myers CR: The developmental dynamics of the maize leaf transcriptome. Nat Genet. 2010, 42 (12): 1060-1067. 10.1038/ng.703.View ArticlePubMed
- Marquez Y, Brown JW, Simpson C, Barta A, Kalyna M: Transcriptome survey reveals increased complexity of the alternative splicing landscape in Arabidopsis. Genome Res. 2012, 22 (6): 1184-1195. 10.1101/gr.134106.111.PubMed CentralView ArticlePubMed
- Filichkin SA, Priest HD, Givan SA, Shen R, Bryant DW, Fox SE, Wong WK, Mockler TC: Genome-wide mapping of alternative splicing in Arabidopsis thaliana. Genome Res. 2010, 20 (1): 45-58. 10.1101/gr.093302.109.PubMed CentralView ArticlePubMed
- Pan Q, Shai O, Lee LJ, Frey BJ, Blencowe BJ: Deep surveying of alternative splicing complexity in the human transcriptome by high-throughput sequencing. Nat Genet. 2008, 40 (12): 1413-1415. 10.1038/ng.259.View ArticlePubMed
- Wang ET, Sandberg R, Luo S, Khrebtukova I, Zhang L, Mayr C, Kingsmore SF, Schroth GP, Burge CB: Alternative isoform regulation in human tissue transcriptomes. Nature. 2008, 456 (7221): 470-476. 10.1038/nature07509.PubMed CentralView ArticlePubMed
- Wang BB, Brendel V: Genomewide comparative analysis of alternative splicing in plants. Proc Natl Acad Sci USA. 2006, 103 (18): 7175-7180. 10.1073/pnas.0602039103.PubMed CentralView ArticlePubMed
- Baek JM, Han P, Iandolino A, Cook DR: Characterization and comparison of intron structure and alternative splicing between Medicago truncatula, Populus trichocarpa Arabidopsis and rice. Plant Mol Biol. 2008, 67 (5): 499-510. 10.1007/s11103-008-9334-4.View ArticlePubMed
- Labadorf A, Link A, Rogers MF, Thomas J, Reddy AS, Ben-Hur A: Genome-wide analysis of alternative splicing in Chlamydomonas reinhardtii. BMC Genomics. 2010, 11: 114-10.1186/1471-2164-11-114.PubMed CentralView ArticlePubMed
- Walters B, Lum G, Sablok G, Min XJ: Genome-wide landscape of alternative splicing events in Brachypodium distachyon. DNA Res. 2013, 20 (2): 163-171. 10.1093/dnares/dss041.PubMed CentralView ArticlePubMed
- Kwan T, Benovoy D, Dias C, Gurd S, Provencher C, Beaulieu P, Hudson TJ, Sladek R, Majewski J: Genome-wide analysis of transcript isoform variation in humans. Nat Genet. 2008, 40 (2): 225-231. 10.1038/ng.2007.57.View ArticlePubMed
- Nakai K, Sakamoto H: Construction of a novel database containing aberrant splicing mutations of mammalian genes. Gene. 1994, 141 (2): 171-177. 10.1016/0378-1119(94)90567-3.View ArticlePubMed
- Yang Z: Small GTPases: versatile signaling switches in plants. Plant Cell. 2002, 14 (Suppl): S375-S388.PubMed CentralPubMed
- Filichkin SA, Mockler TC: Unproductive alternative splicing and nonsense mRNAs: a widespread phenomenon among plant circadian clock genes. Biol Direct. 2012, 7: 20-10.1186/1745-6150-7-20.PubMed CentralView ArticlePubMed
- Li J, Li X, Guo L, Lu F, Feng X, He K, Wei L, Chen Z, Qu LJ, Gu H: A subgroup of MYB transcription factor genes undergoes highly conserved alternative splicing in Arabidopsis and rice. J Exp Bot. 2006, 57 (6): 1263-1273. 10.1093/jxb/erj094.View ArticlePubMed
- Severing EI, van Dijk AD, Morabito G, Busscher-Lange J, Immink RG, van Ham RC: Predicting the impact of alternative splicing on plant MADS domain protein function. PLoS One. 2012, 7 (1): e30524-10.1371/journal.pone.0030524.PubMed CentralView ArticlePubMed
- Craig HL, Wirtz J, Bamps S, Dolphin CT, Hope IA: The significance of alternative transcripts for Caenorhabditis elegans transcription factor genes, based on expression pattern analysis. BMC Genomics. 2013, 14 (1): 249-10.1186/1471-2164-14-249.PubMed CentralView ArticlePubMed
- Li Q, Lin YC, Sun YH, Song J, Chen H, Zhang XH, Sederoff RR, Chiang VL: Splice variant of the SND1 transcription factor is a dominant negative of SND1 members and their regulation in Populus trichocarpa. Proc Natl Acad Sci USA. 2012, 109 (36): 14699-14704. 10.1073/pnas.1212977109.PubMed CentralView ArticlePubMed
- Kalsotra A, Cooper TA: Functional consequences of developmentally regulated alternative splicing. Nat Rev Genet. 2011, 12 (10): 715-729. 10.1038/nrg3052.PubMed CentralView ArticlePubMed
- Blanc G, Wolfe KH: Functional divergence of duplicated genes formed by polyploidy during Arabidopsis evolution. Plant Cell. 2004, 16 (7): 1679-1691. 10.1105/tpc.021410.PubMed CentralView ArticlePubMed
- Freeling M: Bias in plant gene content following different sorts of duplication: tandem, whole-genome, segmental, or by transposition. Annu Rev Plant Biol. 2009, 60: 433-453. 10.1146/annurev.arplant.043008.092122.View ArticlePubMed
- Maere S, De Bodt S, Raes J, Casneuf T, Van Montagu M, Kuiper M, Van de Peer Y: Modeling gene and genome duplications in eukaryotes. Proc Natl Acad Sci USA. 2005, 102 (15): 5454-5459. 10.1073/pnas.0501102102.PubMed CentralView ArticlePubMed
- Paterson AH, Chapman BA, Kissinger JC, Bowers JE, Feltus FA, Estill JC: Many gene and domain families have convergent fates following independent whole-genome duplication events in Arabidopsis, Oryza Saccharomyces and Tetraodon. Trends Genet. 2006, 22 (11): 597-602. 10.1016/j.tig.2006.09.003.View ArticlePubMed
- Birchler JA, Riddle NC, Auger DL, Veitia RA: Dosage balance in gene regulation: biological implications. Trends Genet. 2005, 21 (4): 219-226. 10.1016/j.tig.2005.02.010.View ArticlePubMed
- Thomas BC, Pedersen B, Freeling M: Following tetraploidy in an Arabidopsis ancestor, genes were removed preferentially from one homeolog leaving clusters enriched in dose-sensitive genes. Genome Res. 2006, 16 (7): 934-946. 10.1101/gr.4708406.PubMed CentralView ArticlePubMed
- Nilsen TW, Graveley BR: Expansion of the eukaryotic proteome by alternative splicing. Nature. 2010, 463 (7280): 457-463. 10.1038/nature08909.PubMed CentralView ArticlePubMed
- Zhou R, Moshgabadi N, Adams KL: Extensive changes to alternative splicing patterns following allopolyploidy in natural and resynthesized polyploids. Proc Natl Acad Sci USA. 2011, 108 (38): 16122-16127. 10.1073/pnas.1109551108.PubMed CentralView ArticlePubMed
- Zhang PG, Huang SZ, Pin AL, Adams KL: Extensive divergence in alternative splicing patterns after gene and genome duplication during the evolutionary history of Arabidopsis. Mol Biol Evol. 2010, 27 (7): 1686-1697. 10.1093/molbev/msq054.View ArticlePubMed
- Patel RK, Jain M: NGS QC Toolkit: a toolkit for quality control of next generation sequencing data. PLoS One. 2012, 7 (2): e30619-10.1371/journal.pone.0030619.PubMed CentralView ArticlePubMed
- Benjamini Y, Hochberg Y: Controlling the false discovery rate: a practical and powerful approach to multiple testing. J Roy Statist Soc Ser. 1995, 57 (1): 289-300.
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 (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.