Fast gapped-read alignment with Bowtie 2

Aligning sequencing reads to a reference genome is the first step in many comparative genomics pipelines, including pipelines for variant calling¹, isoform quantitation² and differential gene expression³. In many cases, the alignment step is the slowest. This is because for each read the aligner must solve a difficult computational problem: determining the read's likely point of origin with respect to a reference genome⁴.

Many aligners use a genome index to rapidly narrow the list of candidate alignment locations. The full-text minute index⁵ is a fast and memory-efficient index that has been used in recent aligners^6–10. Index-assisted aligners work by searching for all ways of mutating the read string into a string that occurs in the reference, subject to an alignment policy limiting the number of differences. Although this search space is large, many portions of it can be skipped (‘pruned’) without loss of sensitivity. In practice, pruning strategies such as double indexing⁶ and bidirectional Burrows-Wheeler transform (BWT)⁷ facilitate very efficient ungapped alignment of short reads.

Index-aided alignment can be quite inefficient, however, when alignments are permitted to contain gaps. Alignment gaps can result either from sequencing errors or from true insertions and deletions. Ungapped aligners such as Bowtie will usually fail to align reads spanning gaps and will therefore miss evidence for these events. Gaps greatly increase the size of the search space and reduce the effectiveness of pruning, thereby substantially slowing aligners built solely on index-assisted alignment. Bowtie 2 extends the full-text minute index–based approach of Bowtie to permit gapped alignment by dividing the algorithm broadly into two stages: an initial, ungapped seed-finding stage that benefits from the speed and memory efficiency of the full-text minute index and a gapped extension stage that uses dynamic programming and benefits from the efficiency of single-instruction multiple-data (SIMD) parallel processing available on modern processors. The combination of full-text minute index–assisted seed alignment and SIMD-accelerated dynamic programming achieves an effective combination of speed, sensitivity and accuracy across a range of read lengths and sequencing technologies.

For each read, Bowtie 2 proceeds in four steps (Supplementary Note and Supplementary Fig. 1). In step 1, Bowtie 2 extracts ‘seed’ substrings from the read and its reverse complement. In step 2, the extracted substrings are aligned to the reference in an ungapped fashion assisted by the full-text minute index. In step 3, seed alignments are prioritized, and their positions in the reference genome are calculated from the index. In step 4, seeds are extended into full alignments by performing SIMD-accelerated dynamic programming.

To assess how Bowtie 2 performs on real data, we compared Bowtie 2 to three other full-text minute index–based read aligners: Burrows-Wheeler Aligner (BWA)⁸, BWA's Smith-Waterman alignment (BWA-SW)⁹ and short oligonucleotide alignment program 2 (SOAP2)¹⁰ as well as to Bowtie⁶. In all experiments, the reference we used was the GRCh37 major build of the human genome, including sex chromosomes, mitochondrial genome and ‘non-chromosomal’ sequences. We obtained 100-by-100 nucleotide (nt) paired-end HiSeq (2000) reads from a human resequencing study¹¹ and extracted a random subset of 2 million pairs.

We first used BWA, SOAP2, Bowtie 2 and Bowtie to align one end (labeled ‘1’) from the subset in an unpaired fashion. To illustrate parameter tradeoffs, we ran three of the tools with a wide variety of parameter settings (Fig. 1a and Supplementary Table 1). Note that SOAP2 and Bowtie do not permit gapped alignment of unpaired reads. The Bowtie 2 default mode is faster than all BWA modes we tried and more than 2.5 times faster than the BWA default mode. All Bowtie 2 modes aligned a greater number of reads than either BWA (Supplementary Table 2) or SOAP2. The peak memory footprint of Bowtie 2 (3.24 gigabytes) was between that of BWA (2.39 gigabytes) and SOAP2 (5.34 gigabytes).

We then aligned reads in a paired-end fashion, using a variety of alignment parameters (Fig. 1b and Supplementary Table 1). Bowtie is at a disadvantage in this scenario because it only searches for ungapped, concordant paired-end alignments. We found that the Bowtie 2 default mode was faster than all BWA modes we tried and more than 3 times faster than the BWA default mode. All Bowtie 2 modes aligned a greater number of reads than either BWA (Supplementary Table 2) or SOAP2. The peak memory footprint of Bowtie 2 (3.26 gigabytes) was similar to that of BWA (3.20 gigabytes) and smaller than that of SOAP2 (5.34 gigabytes).

To assess Bowtie 2 performance on longer reads, we obtained Roche 454 reads from the 1000 Genomes Project Pilot¹² and Ion Torrent reads from the G. Moore genome resequencing project¹³. We extracted a random subset of 1 million reads from each and aligned them with BWA-SW and Bowtie 2. We did not align with Bowtie, BWA or SOAP2 because those tools are designed for shorter reads. We configured Bowtie 2 to perform local alignment similar to BWA-SW, and we ran both tools with various parameter settings. For both the 454 and Ion Torrent data (Fig. 1c,d and Supplementary Table 1), the Bowtie 2 default local-alignment mode was faster and aligned more reads (Supplementary Table 2) than any of the BWA-SW modes, with a smaller peak memory footprint (3.39 gigabytes for Bowtie 2 and 3.66 gigabytes for BWA-SW).

To assess the accuracy and sensitivity of Bowtie 2, we used simulated reads for which we knew the correct alignment. Using Mason (http://www.seqan.de/projects/mason.html), we simulated sets of 100,000 Illumina-like single reads 100 nt long and 150 nt long from the human genome or the same number of paired-end reads, and ran Bowtie 2, BWA and SOAP2 on each dataset. We also ran Bowtie on the 100 nt and 100 nt × 100 nt datasets. For each aligner and each dataset, we recorded the number of correct and incorrect alignments stratified by mapping quality, defined as –10 log10(p), where p is the aligner's estimate of the probability that the read was aligned incorrectly. We then calculated the cumulative number of correct and incorrect alignments from high to low mapping quality. We considered an alignment correct only if the leftmost position was within 50 nt of the position assigned by the simulator on the same strand.

We plotted cumulative correct alignments against cumulative incorrect alignments for each dataset and aligner (Fig. 2 and Supplementary Table 3). In all cases, Bowtie 2 and BWA reported more correct alignments than SOAP2 and Bowtie. For the unpaired reads, the plots indicate that Bowtie 2 gave more correct and fewer incorrect alignments than BWA over a range of mapping quality cutoffs. For paired-end reads, the difference was smaller. Note that in its paired-end mode, BWA performed local alignment to recover one of the two ends of the paired-end read in some situations, which Bowtie 2 did not. In the 150-nt dataset comparison, for instance, BWA trimmed 2,991 reads in this way. In Supplementary Results and Supplementary Figure 2 we show results for additional read lengths.

We also used the Mason simulator to generate two collections of 100,000 454-like reads with average lengths of 250 nt and 400 nt. We ran Bowtie 2 and BWA-SW on these datasets (Fig. 2 and Supplementary Table 3). Bowtie 2 generally outperformed BWA-SW, especially for the 250-nt reads. BWA-SW also trimmed more reads; for the 250-nt data, for example, it trimmed 53,486 reads compared to 51,051 reads by Bowtie 2.

Full-text minute index–assisted search is an increasingly common approach for aligning sequencing reads. Extending this method to perform sensitive gapped alignment without incurring serious computational penalties is a major technical challenge. We found that Bowtie 2, a method that combines the advantages of the full-text minute index and SIMD dynamic programming, achieved very fast and memory-efficient gapped alignment of sequencing reads. Bowtie 2 improved on the previous Bowtie method in terms of speed and fraction of reads aligned (Supplementary Results, Supplementary Figs. 3, 4 and Supplementary Tables 4, 5) and was substantially faster than non–full-text minute index–based approaches while aligning a comparable fraction of reads (Supplementary Results and Supplementary Table 6). Robustness to edits, especially gaps, will continue to be a crucial concern as errors typically manifest as gaps in emerging single-molecule sequencing technologies. The speed, sensitivity, accuracy, quality-value awareness and ability to align in both local and end-to-end modes make Bowtie 2 particularly apt for current and future sequencing workloads. Bowtie 2 is free, open-source software available as Supplementary Software and at http://bowtie-bio.sourceforge.net/bowtie2/index.shtml.