tophat bug

[MeixiaZhao]

Hello,

It seemed I found a bug in the previous tophat versions which use bowtie as the aligment tool. I tested in this way:

Test1: tophat using bowtie

Option1: tophat -r 0 -i 30 -I 100000 -p 8 -g 1 -o g1 reference R1_101.fastq R2_101.fastq;

Option2: tophat -r 0 -i 30 -I 100000 -p 8 -g 10 -o g1 reference R1_101.fastq R2_101.fastq;

I caculated the percentages of unique reads for both Option1 and Option2. Actually, all the reads in Option1 are unique, the percentage is 60.8%. For Option2, I only caculated the unique ones, whose percentage is 82.3% (I used samtools to transform the .bam file into .sam file and caculated unique reads). Firstly, this made me so confusing, why this made so large differences. So later on I compared the two .sam files and I found even some unique reads really exist, when you use -g 1, you will loss some data. I searched some web comments, they said that bowtie has no gap-alignment and is not good for the reads longer than 50bp.

Test2: tophat using bowtie2

Option3: tophat -r 0 -i 30 -I 100000 -p 8 -g 1 -o g1_bowtie2 reference R1_101.fastq R2_101.fastq;

Option4: tophat -r 0 -i 30 -I 100000 -p 8 -g 10 -o g1_bowtie2 reference R1_101.fastq R2_101.fastq;

For Option3 the unique reads percentage is 82.9%, Option4 the unique reads percentage is 83.8%, only 0.9% difference. It seems bowtie2 did much better job than bowtie.

Did anyone meet this kind of situation? If you have, please just let me know. I'm just worring about I did something wrong. I appreciated your kindness.

[sdriscoll]

Option1 and Option2 seem to be the same command. Is Option2 supposed to read -g 10?

[MeixiaZhao]

Yes, I made a mistake.

[sdriscoll]

i think what you are seeing here is fine. I'd expect bowtie2 to be better at aligning longer reads than bowtie[1].

with respect to the differences you see in option 1 verses option 2 the reason is the "behind the scenes" options passed to bowtie. When you run -g 1 in tophat bowtie receives the options -m 1 -k 1. the -m option in bowtie says "don't report alignments for any reads that can align more than N times" and the -k option says "report up to N alignments per read". You could think of -k as redundant when -m is specified in this way. When you use -g 10 then bowtie gets -m 10 -k 10 so it allows reads that can align up to 10 times to be reported.

I guess the strange thing is that the -g 10 option produced more unique alignments than -g 1. This should not be the case. It would make sense if bowtie was using the --best --strata options but it's not. In fact regardless of the -g setting you should see the exact same number of unique alignments.

Maybe what's happening is a result of the paired reads. Maybe some reads that could align in multiple locations only had one unique, valid, paired alignment. That could explain it.

[MeixiaZhao]

Actually, I just want to use unique reads to do later analysis. I chose Option3, as it didn't lose ao much data. Am I right?

[sdriscoll]

Yes, go with whatever produces the most unique alignments do you are using as much of your data as possible.