Efficient RNA isoform identification and quantification from RNA-Seq data with network flows

Abstract

Motivation: Several state-of-the-art methods for isoform identification and quantification are based on [Formula: see text]-regularized regression, such as the Lasso. However, explicitly listing the-possibly exponentially-large set of candidate transcripts is intractable for genes with many exons. For this reason, existing approaches using the [Formula: see text]-penalty are either restricted to genes with few exons or only run the regression algorithm on a small set of preselected isoforms.

Results: We introduce a new technique called FlipFlop, which can efficiently tackle the sparse estimation problem on the full set of candidate isoforms by using network flow optimization. Our technique removes the need of a preselection step, leading to better isoform identification while keeping a low computational cost. Experiments with synthetic and real RNA-Seq data confirm that our approach is more accurate than alternative methods and one of the fastest available.

Availability and implementation: Source code is freely available as an R package from the Bioconductor Web site (http://www.bioconductor.org/), and more information is available at http://cbio.ensmp.fr/flipflop.

Supplementary information: Supplementary data are available at Bioinformatics online.

Fig. 1.

Illustration of the graph construction for a gene with 5 exons. The original splicing graph is represented in (a). The 5 exons are represented as vertices and an arrow between two vertices indicates a junction. The nodes of graph $G\prime$ in (b) and (c) are bins with positive effective length denoted by gray square, as well as source s and sink t represented as circles. $G\prime$ in (b) is the resulting graph when all exons are bigger than the read length. In that case, each bin either corresponds to a unique exon, or to a junction between two exons. $G\prime$ in (c) is the resulting graph when the length of exon 3 is smaller than the read length. Some bins involve then more than two exons, here bins $(2\text{-}3\text{-}4)$ and $(2\text{-}3\text{-}5)$. The source links all possible starting bins, and conversely all possible stopping bins are linked to the sink. There is a one-to-one correspondence between (s, t)-paths in $G\prime$ (paths starting at s and ending at t) and isoform candidates. For example, the path $(s,1,1\text{-}4,4,4\text{-}5,5,t)$ corresponds to isoform 1-4-5

Fig. 2.

Flow interpretation of isoforms using the same graph as in Figure 1b. For the sake of clarity, some edges connecting s and t to internal nodes are not represented, and the length of the different bins is assumed to be equal. In (a), one unit of flow is carried along the path in red, corresponding to an isoform with abundance 1. In (b), another isoform with abundance 3 is added, yielding additional read counts at every node

Fig. 3.

Precision and recall of compared methods on simulated reads from the UCSC annotated human transcripts. (a) Single-end reads with different lengths (100, 200, 300 bp). (b) Paired-end reads with different lengths (100, 125, 150, 175 bp). (c) Single-end reads with a fixed 150 bp length and an increasing amout of material (1, 5, 10 million)

Fig. 4.

Average CPU times in milliseconds (logarithmic scale) for the compared methods with process a gene from human simulated 100 bp single-end reads

Fig. 5.

Precision and recall of compared methods on human embryonic stem cells data