RNA-Seq

RNA-seq (RNA sequencing), also called whole transcriptome shotgun sequencing (WTSS), uses next-generation sequencing (NGS) to reveal the presence and quantity of RNA in a biological sample at a given moment in time.

RNA-Seq is used to analyze the continually changing cellular transcriptome. Specifically, RNA-Seq facilitates the ability to look at alternative gene spliced transcripts, post-transcriptional modifications, gene fusion, mutations/SNPs and changes in gene expression. In addition to mRNA transcripts, RNA-Seq can look at different populations of RNA to include total RNA, small RNA, such as miRNA, tRNA, and ribosomal profiling. RNA-Seq can also be used to determine exon/intron boundaries and verify or amend previously annotated 5’ and 3’ gene boundaries.

Prior to RNA-Seq, gene expression studies were done with hybridization-based microarrays. Issues with microarrays include cross-hybridization artifacts, poor quantification of lowly and highly expressed genes, and the knowledge of the sequence. Because of these technical issues, transcriptomics transitioned to sequencing-based methods. These progressed from Sanger sequencing of Expressed Sequence Tag libraries, to chemical tag-based methods (e.g., serial analysis of gene expression), and finally to the current technology, NGS of cDNA (notably RNA-Seq).

Methods

File:RNASeqPics1.jpg

File:RNASeqPics2.jpg

RNA 'Poly(A)' library

See also: Polyadenylation

Creation of a sequence library can change from platform to platform in high throughput sequencing, where each has several kits designed to build different types of libraries and adapting the resulting sequences to the specific requirements of their instruments. However, due to the nature of the template being analyzed, there are commonalities within each technology. Frequently, in mRNA analysis the 3' polyadenylated (poly(A)) tail is targeted in order to ensure that coding RNA is separated from noncoding RNA. This can be accomplished simply with poly (T) oligos covalently attached to a given substrate. Presently many studies utilize magnetic beads for this step.

Studies including portions of the transcriptome outside poly(A) RNAs have shown that when using poly(T) magnetic beads, the flow-through RNA (non-poly(A) RNA) can yield important noncoding RNA gene discovery which would have otherwise gone unnoticed. Also, since ribosomal RNA represents over 90% of the RNA within a given cell, studies have shown that its removal via probe hybridization increases the capacity to retrieve data from the remaining portion of the transcriptome.

The next step is reverse transcription. Due to the 5' bias of randomly primed-reverse transcription as well as secondary structures influencing primer binding sites, hydrolysis of RNA into 200-300 nucleotides prior to reverse transcription reduces both problems simultaneously. However, there are trade-offs with this method where although the overall body of the transcripts are efficiently converted to DNA, the 5' and 3' ends are less so. Depending on the aim of the study, researchers may choose to apply or ignore this step.

Once the cDNA is synthesized it can be further fragmented to reach the desired fragment length of the sequencing system.

Small RNA/non-coding RNA sequencing

When sequencing RNA other than mRNA, the library preparation is modified. The cellular RNA is selected based on the desired size range. For small RNA targets, such as miRNA, the RNA is isolated through size selection. This can be performed with a size exclusion gel, through size selection magnetic beads, or with a commercially developed kit. Once isolated, linkers are added to the 3’ and 5’ end then purified. The final step is cDNA generation through reverse transcription.

File:RNA-Seq-alignment.png

RNA-seq mapping of short reads in exon-exon junctions.

Direct RNA sequencing

As converting RNA into cDNA using reverse transcriptase has been shown to introduce biases and artifacts that may interfere with both the proper characterization and quantification of transcripts, single molecule Direct RNA Sequencing (DRSTM) technology was under development by Helicos (now bankrupt). DRSTM sequences RNA molecules directly in a massively-parallel manner without RNA conversion to cDNA or other biasing sample manipulations such as ligation and amplification.

Transcriptome assembly

See also: Sequence alignment software § Short-Read Sequence Alignment

Two different assembly methods are used for producing a transcriptome from raw sequence reads: de-novo and genome-guided.

The first approach does not rely on the presence of a reference genome in order to reconstruct the nucleotide sequence. Due to the small size of the short reads, de novo assembly may be difficult, though some software does exist ( Velvet (algorithm), Oases, and Trinity to mention a few), as there cannot be large overlaps between each read needed to easily reconstruct the original sequences. The deep coverage also makes the computing power to track all the possible alignments prohibitive. This deficit can be improved using longer sequences obtained from the same sample using other techniques such as Sanger sequencing, and using larger reads as a "skeleton" or a "template" to help assemble reads in difficult regions (e.g. regions with repetitive sequences).

An “easier” and relatively computationally cheaper approach is that of aligning the millions of reads to a "reference genome". There are many tools available for aligning genomic reads to a reference genome (sequence alignment tools), however, special attention is needed when aligning a transcriptome to a genome, mainly when dealing with genes having intronic regions. Several software packages exist for short read alignment, and recently specialized algorithms for transcriptome alignment have been developed, e.g. Bowtie for RNA-seq short read alignment, TopHat for aligning reads to a reference genome to discover splice sites, Cufflinks to assemble the transcripts and compare/merge them with others, or FANSe. These tools can also be combined to form a comprehensive system.

Although numerous solutions to the assembly quest have been proposed, there is still lots of room for improvement given the resulting variability of the approaches. A group from the Center for Computational Biology at the East China Normal University in Shanghai compared different de novo and genome-guided approaches for RNA-Seq assembly. They noted that, although most of the problems can be solved using graph theory approaches, there is still a consistent level of variability in all of them. Some algorithms outperformed the common standards for some species while still struggling for others. The authors suggest that the “most reliable” assembly could be then obtained by combining different approaches. Interestingly, these results are consistent with NGS-genome data obtained in a recent contest called Assemblathon where 21 contestants analyzed sequencing data from three different vertebrates (fish, snake and bird) and handed in a total of 43 assemblies. Using a metric made of 100 different measures for each assembly, the reviewers concluded that 1) assembly quality can vary a lot depending on which metric is used and 2) assemblies that scored well in one species did not really perform well in the other species.

As discussed above, sequence libraries are created by extracting mRNA using its poly(A) tail, which is added to the mRNA molecule post-transcriptionally and thus splicing has taken place. Therefore, the created library and the short reads obtained cannot come from intronic sequences, so library reads spanning the junction of two or more exons will not align to the genome.

A possible method to work around this is to try to align the unaligned short reads using a proxy genome generated with known exonic sequences. This need not cover whole exons, only enough so that the short reads can match on both sides of the exon-exon junction with minimum overlap. Some experimental protocols allow the production of strand specific reads.

Experimental considerations

The information gathered when sequencing a sample's transcriptome in this way has many of the same limitations and advantages as other RNA expression analysis pipelines. The main pros and cons of this approach can be summarized as:

a) Tissue specificity: Gene expression is not uniform throughout an organism's cells, it is strongly dependent on the tissue type being measured; RNA-Seq, as any other sequencing technology that analyzes homogeneous samples, can provide a complete snapshot of all the transcripts being available at that precise moment in the cell. This approach is unlikely to be biased like an oligonucleotide microarray approach that instead analyzes a selected number of previously defined transcripts.

b) Time dependent: During a cell's lifetime and context, its gene expression levels change. As previously mentioned any single sequencing experiment will offer information regarding one point in time. Time course experiments are so far the only solution that would allow a complete overview of the circadian transcriptome so that researchers could obtain a precise description of the physiological changes happening over time. However, this approach is unfeasible for patient samples since it is quite improbable that biopsies will be collected serially in short time intervals. A possible work-around could be the use of urine, blood or saliva samples that won’t require any invasive procedure.

c) Coverage: coverage/depth can affect the mutations seen. Given that everything is expression-centric, an allele might not be detected, either because it is not in the genome, or because it is not being expressed. At the same time, RNA-seq can yield additional information rather than just the existence of a heterozygous gene as it can also help in estimating the expression of each allele. In association studies, genotypes are associated to disease and expression levels can also be associated with disease. Using RNA-seq, we can measure the relationship between these two associated variables, that is, in what relation are each of the alleles being expressed.

The depth of sequencing required for specific applications can be extrapolated from a pilot experiment.

d) Subjectivity of the analysis: As described above, numerous attempts have been taken to uniformly analyze the data. However, the results can vary due to the multitude of algorithms and pipelines available. Most of the approaches are correct, but have to be tailored to the needs of the investigators in order to better capture the desired effect. This variability in methods, although in smaller scale, is still present in other RNA profiling approaches where reagents, personnel and techniques can lead to similar, although statistically different, results. Because of this, care must be taken when drawing conclusions from the sequencing experiment, as some information gathered might not be representative of the individual.

e) Data management: The main issue with NGS data is the volume of data produced. Microarray data occupy up to one thousand times less disk space than NGS data therefore requiring smaller storage units. The high capacity storage units required by RNA-Seq data are, however, directly proportional to the volume of information that goes with it. The payoff of “more complete” big scale datasets have to be evaluated prior to starting the experiment.

f) Downstream interpretation of the data: Different layers of interpretations have to be considered when analyzing RNA-Seq data. Biological, clinical and regulatory functions of the results are what allow clinicians and investigators to draw meaningful conclusions (i.e. the sequence of an RNA molecule presents, although identified with different read depths, might not perfectly mirror the initial DNA sequence). An example of this would be during SNV discovery as the mutations discovered are more precisely the mutations being expressed. Observing a homozygote location to a non-reference allele in an organism does not necessarily mean that this is the individual's genotype, it could just mean that the gene copy with the reference allele is not being expressed in that tissue and/or at the time snapshot the sample was acquired.

Analysis

See also: List of RNA-Seq bioinformatics tools

Gene expression

The characterization of gene expression in cells via measurement of mRNA levels has long been of interest to researchers, both in terms of which genes are expressed in what tissues, and at what levels. Even though it has been shown that due to other post transcriptional gene regulation events (such as RNA interference) there is not necessarily always a strong correlation between the abundance of mRNA and the related proteins, measuring mRNA concentration levels is still a useful tool in determining how the transcriptional machinery of the cell is affected in the presence of external signals (e.g. drug treatment), or how cells differ between a healthy state and a diseased state.

Expression can be deduced via RNA-seq to the extent at which a sequence is retrieved. Transcriptome studies in yeast show that in this experimental setting, a fourfold coverage is required for amplicons to be classified and characterized as an expressed gene. When the transcriptome is fragmented prior to cDNA synthesis, the number of reads corresponding to the particular exon normalized by its length in vivo yields gene expression levels which correlate with those obtained through qPCR. This is frequently further normalized by the total number of mapped reads so that expression levels are expressed as Fragments Per Kilobase of transcript per Million mapped reads (FPKM).

The only way to be absolutely sure of the individual's mutations is to compare the transcriptome sequences to the germline DNA sequence. This enables the distinction of homozygous genes versus skewed expression of one of the alleles and it can also provide information about genes that were not expressed in the transcriptomic experiment. An R-based statistical package known as CummeRbund can be used to generate expression comparison charts for visual analysis.

Differential expression and absolute quantification of transcripts

RNA-Seq is generally used to compare gene expression between conditions, such as a drug treatment vs non-treated, and find out which genes are up- or down regulated in each condition. In principle RNA-Seq will allow to account for all the transcripts in the cell for each condition. Differently expressed genes can be identified by using tools that count the sequencing reads per gene and compare them between samples. The most commonly used tools for this type of analysis are DESeq and edgeR, packages from Bioconductor. Both these tools use a model based on the negative binomial distribution.

It is not possible to do absolute quantification using the common RNA-Seq pipeline, because it only provides RNA levels relative to all transcripts. If the total amount of RNA in the cell changes between conditions, relative normalization will misrepresent the changes for individual transcripts. Absolute quantification of mRNAs is possible by performing RNA-Seq with added spike ins, samples of RNA at known concentrations. After sequencing, the read count of the spike ins sequences is used to determine the direct correspondence between read count and biological fragments. In developmental studies, this technique has been used in Xenopus tropicalis embryos at a high temporal resolution, to determine transcription kinetics.

Single nucleotide variation discovery

See also: single nucleotide polymorphism

Transcriptome single nucleotide variation has been analyzed in maize on the Roche 454 sequencing platform. Directly from the transcriptome analysis, around 7000 single nucleotide polymorphisms (SNPs) were recognized. Following Sanger sequence validation, the researchers were able to conservatively obtain almost 5000 valid SNPs covering more than 2400 maize genes. RNA-seq is limited to transcribed regions however, since it will only discover sequence variations in exon regions. This misses many subtle but important intron alleles that affect disease such as transcription regulators, leaving analysis to only large effectors. While some correlation exists between exon to intron variation, only whole genome sequencing would be able to capture the source of all relevant SNPs.

Post-transcriptional SNVs

Having the matching genomic and transcriptomic sequences of an individual can also help in detecting post-transcriptional edits, where, if the individual is homozygous for a gene, but the gene's transcript has a different allele, then a post-transcriptional modification event is determined.

mRNA centric single nucleotide variants (SNVs) are generally not considered as a representative source of functional variation in cells, mainly due to the fact that these mutations disappear with the mRNA molecule, however the fact that efficient DNA correction mechanisms do not apply to RNA molecules can cause them to appear more often. This has been proposed as the source of certain prion diseases, also known as TSE or transmissible spongiform encephalopathies.

File:RNA-Seq-fusion-gene.png

RNA-seq mapping of short reads over exon-exon junctions, depending on where each end maps to, it could be defined a Trans or a Cis event.

Fusion gene detection

See also: Fusion gene

Caused by different structural modifications in the genome, fusion genes have gained attention because of their relationship with cancer. The ability of RNA-seq to analyze a sample's whole transcriptome in an unbiased fashion makes it an attractive tool to find these kinds of common events in cancer.

The idea follows from the process of aligning the short transcriptomic reads to a reference genome. Most of the short reads will fall within one complete exon, and a smaller but still large set would be expected to map to known exon-exon junctions. The remaining unmapped short reads would then be further analyzed to determine whether they match an exon-exon junction where the exons come from different genes. This would be evidence of a possible fusion event, however, because of the length of the reads, this could prove to be very noisy. An alternative approach is to use pair-end reads, when a potentially large number of paired reads would map each end to a different exon, giving better coverage of these events (see figure). Nonetheless, the end result consists of multiple and potentially novel combinations of genes providing an ideal starting point for further validation.

Coexpression networks

Coexpression networks are data-derived representations of genes behaving in a similar way across tissues and experimental conditions. Their main purpose lies in hypothesis generation and guilt-by-association approaches for inferring functions of previously unknown genes. RNASeq data has been recently used to infer genes involved in specific pathways based on Pearson correlation, both in plants and mammals. The main advantage of RNASeq data in this kind of analysis over the microarray platforms is the capability to cover the entire transcriptome, therefore allowing the possibility to unravel more complete representations of the gene regulatory networks. Differential regulation of the splice isoforms of the same gene can be detected and used to predict and their biological functions. Weighted gene co-expression network analysis has been successfully used to identify co-expression modules and intramodular hub genes based on RNA seq data. Co-expression modules may corresponds to cell types or pathways. Highly connected intramodular hubs can be interpreted as representatives of their respective module. Variance-Stabilizing Transformation approaches for estimating correlation coefficients based on RNA seq data have been proposed.

Application to genomic medicine

History

The past five years have seen a flourishing of NGS-based methods for genome analysis leading to the discovery of a number of new mutations and fusion transcripts in cancer. RNA-Seq data could help researchers interpreting the “personalized transcriptome” so that it will help understanding the transcriptomic changes happening therefore, ideally, identifying gene drivers for a disease. The feasibility of this approach is however dictated by the costs in terms of money and time.

A basic search on PubMed reveals that the term RNA Seq, queried as “"RNA Seq" OR "RNA-Seq" OR "RNA sequencing" OR "RNASeq"” in order to capture the most common ways of phrasing it, gives 5.425 hits demonstrating usage statistics of this technology. A few examples will be taken into consideration to explain that RNA-Seq applications to the clinic have the potentials to significantly affect patient’s life and, on the other hand, requires a team of specialists (bioinformaticians, physicians/clinicians, basic researchers, technicians) to fully interpret the huge amount of data generated by this analysis.

As an example of clinical applications, researchers at the Mayo Clinic used an RNA-Seq approach to identify differentially expressed transcripts between oral cancer and normal tissue samples. They also accurately evaluated the allelic imbalance (AI), ratio of the transcripts produced by the single alleles, within a subgroup of genes involved in cell differentiation, adhesion, cell motility and muscle contraction identifying a unique transcriptomic and genomic signature in oral cancer patients. Novel insight on skin cancer (melanoma) also come from RNA-Seq of melanoma patients. This approach led to the identification of eleven novel gene fusion transcripts originated from previously unknown chromosomal rearrangements. Twelve novel chimeric transcripts were also reported, including seven of those that confirmed previously identified data in multiple melanoma samples. Furthermore, this approach is not limited to cancer patients. RNA-Seq has been used to study other important chronic diseases such as Alzheimer (AD) and diabetes. In the former case, Twine and colleagues compared the transcriptome of different lobes of deceased AD’s patient’s brain with the brain of healthy individuals identifying a lower number of splice variants in AD’s patients and differential promoter usage of the APOE-001 and -002 isoforms in AD’s brains. In the latter case, different groups showed the unicity of the beta-cells transcriptome in diabetic patients in terms of transcripts accumulation and differential promoter usage and long non coding RNAs (lncRNAs) signature.

Compared with microarrays, NGS technology has identified novel and low frequency RNAs associated with disease processes. This advantage aids in the diagnosis and possible future treatments of diseases, including cancer. For example, NGS technology identified several previously undocumented differentially-expressed transcripts in rats treated with AFB1, a potent hepatocarcinogen. Nearly 50 new differentially-expressed transcriptions were identified between the controls and AFB1-treated rats. Additionally potential new exons were identified, including some that are responsive to AFB1. The next-generation sequencing pipeline identified more differential gene expressions compared with microarrays, particularly when DESeq software was utilized. Cufflinks identified two novel transcripts that were not previously annotated in the Ensembl database; these transcripts were confirmed using cloning PCR. Numerous other studies have demonstrated NGS's ability to detect aberrant mRNA and small non-coding RNA expression in disease processes above that provided by microarrays. The lower cost and higher throughput offered by NGS confers another advantage to researchers.

The role of small non-coding RNAs in disease processes has also been explored in recent years. For example, Han et al. (2011) examined microRNA expression differences in bladder cancer patients in order to understand how changes and dysregulation in microRNA can influence mRNA expression and function. Several microRNAs were differentially expressed in the bladder cancer patients. Upregulation in the aberrant microRNAs was more common than downregulation in the cancer patients. One of the upregulated microRNAs, hsa-miR-96, has been associated with carcinogenesis, and several of the overexpressed microRNAs have also been observed in other cancers, including ovarian and cervical. Some of the downregulated microRNAs in cancer samples were hypothesized to have inhibitory roles.

ENCODE and TCGA

A lot of emphasis has been given to RNA-Seq data after the Encyclopedia of the regulatory elements (ENCODE) and The Cancer Genome Atlas (TCGA) projects have used this approach to characterize dozens of cell lines and thousands of primary tumor samples, respectively. The former aimed to identify genome-wide regulatory regions in different cohort of cell lines and transcriptomic data are paramount in order to understand the downstream effect of those epigenetic and genetic regulatory layers. The latter project, instead, aimed to collect and analyze thousands of patient’s samples from 30 different tumor types in order to understand the underlying mechanisms of malignant transformation and progression. In this context RNA-Seq data provide a unique snapshot of the transcriptomic status of the disease and look at an unbiased population of transcripts that allows the identification of novel transcripts, fusion transcripts and non-coding RNAs that could be undetected with different technologies.

External links

• RNA-Seq for Everyone: a high-level guide to designing and implementing an RNA-Seq experiment.
• ChIPBase database: provides expression profiles of protein-coding genes and lncRNAs (lincRNAs) from RNA-Seq data across 22 tissues.
• Martin A. Perdacher (September 2011) Next-Generation Sequencing and its Applications in RNA-Seq. Theory part of the Bachelorthesis, Hagenberg.