Mapping short reads

Regular mappers

These mapping tools are useful for reads of DNA origin that should map to a continuous stretch of genomic DNA. Some of these tools can tolerate short indels but they're not designed for reads that span a splice junction

One may choose between bowtie version 1 (faster but ignores indels) and bowtie version 2 (slower but performs gapped alignment (i.e., indels)). For a feature comparision, see ​How is Bowtie 2 different from Bowtie 1?

​bowtie version 1

Sample command:

bsub bowtie  -k 1 -n 2 -l 70 --best --sam --solexa1.3-quals /nfs/genomes/mouse_gp_jul_07_no_random/bowtie/mm9 s_7.txt s7_mm9.k1.n2.l36.best.sam

Parameters included in the sample command:

• -l/--seedlen <int> seed length for -n (default: 28) -- Set to longest possible length of high-quality bases. Use the FastQC output to determine length of high-quality positions.
• -n/--seedmms <int> max mismatches in seed (can be 0-3, default: -n 2)
• -k <int> report up to <int> good alignments per read (default: 1) -- If you want only uniquely mapped reads, however, also use '-m 1' to ignore multi-mapped reads
• --best (in the case of multi-mapped reads, keep only the best hit(s))
• --sam to get SAM output format (which is the best format for downstream analysis)

Choices for fastq encoding (which is listed as "Encoding" in the top "Basic Statistics" table of the FastQC output file). See the ​FASTQ format page for more details.

• --solexa-quals (for input quality scores from Illumina versions 1.2 and earlier)
• --solexa1.3-quals or --phred64-quals (for input quality scores from Illumina versions 1.3-1.7)
• --phred33-quals (default "Sanger format"; for input quality scores from Illumina versions 1.8 and later)

To see other parameters log into tak and type bowtie

​bowtie version 2

Bowtie 2 was designed as an improvement to bowtie 1. See the first ​bowtie2 FAQ for how they differ. Early versions of bowtie 2 had some issues, but these seem to have been fixed. Bowtie 2 uses a different set of genome index files (*.bt2) than bowtie 1 (*.ebwt). Bowtie 2 works with indels.

Sample command:

bsub bowtie2 --phred64 -L 22 -N 1 -x /nfs/genomes/mouse_gp_jul_07_no_random/bowtie/mm9 s_7.txt -S s7_mm9.L22.N1.sam

The parameters included in the sample command:

• -L <int> length of seed substrings; must be >3 and <32 (default=22)
• -N <int> max # mismatches in seed alignment; can be 0 or 1 (default=0)
• -S name of SAM output file

Choices for fastq encoding (which is listed as "Encoding" in the top "Basic Statistics" table of the FastQC output file). See the ​FASTQ format page for more details.

• --solexa-quals (for input quality scores from Illumina versions 1.2 and earlier)
• --phred64 (for input quality scores from Illumina versions 1.3-1.7)
• --phred33 (default "Sanger format"; for input quality scores from Illumina versions 1.8 and later)

bowtie2 can also perform local alignments where the unaligned end(s) of a read are clipped (so, for example, remaining adapter won't prevent alignment) by adding the argument --local.

​bwa - Burrows-Wheeler Alignment Tool

Bwa is a software package containing several related algorithms using the Burrows-Wheeler Transform. It works well even with indels, but not with spliced (RNA) reads.

Sample commands for short (<=100 nt) reads:

# Align single-end reads
bsub "bwa aln /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Reads.fq > Mapped_reads.sai"
bsub "bwa samse /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Mapped_reads.sai Reads.fq > Reads.bwa.sam"
# Align paired-end reads
bsub "bwa aln /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Reads.1.fq > Mapped_reads.1.sai"
bsub "bwa aln /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Reads.2.fq > Mapped_reads.2.sai"
bsub "bwa sampe /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Mapped_reads.1.sai Mapped_reads.2.sai Reads.1.fq Reads.2.fq > Reads.bwa.sam"

Sample commands for long (70 nt - 1M nt) reads:

# Align single-end reads
bsub "bwa mem /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Reads.fq > Mapped_reads.bwa.sam"
# Align paired-end reads
bsub "bwa mem /nfs/genomes/mouse_gp_jul_07_no_random/bwa/mm9 Reads.1.fq Reads.1.fq > Mapped_reads.bwa.sam"

Other tools

Many other regular mapping tools are also available, although they generally require a tool-specific indexed version of the genome.

Splice-aware mappers

These mappers permit the beginning and end of a read to map to (originate from) different places in the genome, which is common for spliced RNA.

​tophat version 1

Running TopHat version 1 requires a change to a user's environment on tak (and only applies to the specific tak session. First run this command:

export PATH="/usr/local/share/tophat1:$PATH"

and then check that your terminal will use the correct TopHat version:

tophat --version

Sample command:

bsub tophat -o s_7_tophat_out --phred64-quals --segment-length 20 -I 200000 -G /nfs/genomes/mouse_gp_jul_07_no_random/gtf/Mus_musculus.NCBIM37.67_noNT.gtf --no-novel-juncs /nfs/genomes/mouse_gp_jul_07_no_random/bowtie/mm9 s_7.txt

The parameters included in the sample command are:

• -o/--output-dir <word> All output files will be created in this directory (default = tophat_out)
• --segment-length <int> Shortest length of a spliced read that can map to one side of the junction. For reads shorter than ~45 nt, set this to half the read length (so set '--segment-length 20' for 40-nt reads). For longer reads, the default length (25) can be used.
• -I <int> Maximum intron length. If your genome has introns that are all shorter (or many that are longer) than the default value (500000), set this to a more appropriate value.
• -G <GTF file> Supply bowtie with a GTF file of transcript models. This can help bowtie identify functions that may otherwise be missed.
• --no-novel-juncs Only look for spliced reads across junctions in the supplied GTF file. Typically not used.

Choices for fastq encoding (which is listed as "Encoding" in the top "Basic Statistics" table of the FastQC output file). See the ​FASTQ format page for more details.

• --solexa-quals (for input quality scores from Illumina versions 1.2 and earlier)
• --solexa1.3-quals or --phred64-quals (for input quality scores from Illumina versions 1.3-1.7)

​tophat version 2

TopHat version 2 uses bowtie2, rather than bowtie, for its mapping. As a result, TopHat 2 uses a different set of genome index files (*.bt2) than TopHat 1 (*.ebwt).

Sample command:

bsub tophat -o s_7_tophat_out --phred64-quals --segment-length 20 -I 200000 -G /nfs/genomes/mouse_gp_jul_07_no_random/gtf/Mus_musculus.NCBIM37.67_noNT.gtf --no-novel-juncs /nfs/genomes/mouse_gp_jul_07_no_random/bowtie/mm9 s_7.txt

The parameters included in the sample command are:

• -o/--output-dir <word> All output files will be created in this directory (default = tophat_out)
• --segment-length <int> Shortest length of a spliced read that can map to one side of the junction. For reads shorter than ~45 nt, set this to half the read length (so set '--segment-length 20' for 40-nt reads). For longer reads, the default length (25) can be used.
• -I <int> Maximum intron length. If your genome has introns that are all shorter (or many that are longer) than the default value (500000), set this to a more appropriate value.
• -G <GTF file> Supply bowtie with a GTF file of transcript models. This can help bowtie identify functions that may otherwise be missed.
• --no-novel-juncs Only look for spliced reads across junctions in the supplied GTF file. Typically not used.

Choices for fastq encoding (which is listed as "Encoding" in the top "Basic Statistics" table of the FastQC output file). See the ​FASTQ format page for more details.

• --solexa-quals (for input quality scores from Illumina versions 1.2 and earlier)
• --solexa1.3-quals or --phred64-quals (for input quality scores from Illumina versions 1.3-1.7)
• For "Sanger / Illumina 1.8" or "Sanger / Illumina 1.9", bowtie can use the default "phred33" encoding