Errors in RNA-Seq quantification affect genes of relevance to human disease

Abstract

Background: RNA-Seq has emerged as the standard for measuring gene expression and is an important technique often used in studies of human disease. Gene expression quantification involves comparison of the sequenced reads to a known genomic or transcriptomic reference. The accuracy of that quantification relies on there being enough unique information in the reads to enable bioinformatics tools to accurately assign the reads to the correct gene.

Results: We apply 12 common methods to estimate gene expression from RNA-Seq data and show that there are hundreds of genes whose expression is underestimated by one or more of those methods. Many of these genes have been implicated in human disease, and we describe their roles. We go on to propose a two-stage analysis of RNA-Seq data in which multi-mapped or ambiguous reads can instead be uniquely assigned to groups of genes. We apply this method to a recently published mouse cancer study, and demonstrate that we can extract relevant biological signal from data that would otherwise have been discarded.

Conclusions: For hundreds of genes in the human genome, RNA-Seq is unable to measure expression accurately. These genes are enriched for gene families, and many of them have been implicated in human disease. We show that it is possible to use data that may otherwise have been discarded to measure group-level expression, and that such data contains biologically relevant information.

Fig. 1

Comparison of methods on global simulated data. a Scatter plots comparing FPKM for each of the 12 methods against the known FPKM from simulated data. The red line indicates the y = x line. b Histograms of read counts for each of the 12 methods. All methods should have a single peak at 1000. c A heatmap of read counts from 843 grossly underestimated genes and 187 grossly overestimated genes. Black and darker colours indicate read counts close to 1000 (accurate); green colours indicate underestimation and red colours overestimation

Fig. 2

General characteristics of problematic genes. Boxplots comparing the length of the shortest exon, the length of the longest exon, the mean exon length, the total number of exons, the transcript length, transcript percentage GC, the number of reads overlapping from the STAR alignment and the number of reads overlapping the TopHat alignment for the 958 problematic genes and the 18,696 other genes

Fig. 3

Comparison of methods on difficult genes. a Scatter plots comparing observed FPKM for each of the 12 methods against the known FPKM from simulated data. The red line indicates the y = x line. b Scatter plots comparing observed read counts for each of the 12 methods against the known read counts from simulated data. The red line indicates the y = x line

Fig. 4

Principle components analysis (PCA) of mouse cancer study. a PCA of tumour (red) and normal (blue) RNA-Seq datasets from each of five cell types. Input data are log(FPKM) values after mapping data using STAR and counting only uniquely mapped reads against known mouse genes (stage 1 analysis) (b) PCA of tumour (red) and normal (blue) RNA-Seq datasets from each of five cell types. Input data are log(FPM) values of reads that cannot be assigned to a single gene but can be uniquely assigned to a multi-map group (MMG). The reads used in (b) are only those reads discarded from (a)

Fig. 5

Heatmap of novel multi-map groups (MMGs). A heatmap of the log FPM (fragments per million) values for 672 differentially expressed MMGs that do not contain any genes present in the list of differentially expressed genes from an analysis of unique counts. The heatmap demonstrates that MMGs which are exclusive of differentially expressed genes from unique counts can be used to separate tumour from normal samples

Fig. 6

Comparison of read counts for (a) ENSMUSG00000024121 and (b) MG4194. Read counts expressed as a percentage of the mapped reads for gene ENSMUSG00000024121, and MG4194, a single-gene MMG that contains only ENSMUSG00000024121. ENSMUSG00000024121 was not found to be differentially expressed by the unique read analysis, but MG4194 was found to be differentially expressed by the MMG analysis. Black bars represent tumour samples, white bars normal samples