Breakdancer Options, Unconsistency Please Help!

Question

[madkitty]

I want to run a CNV call with breakdancer, so I tried two methods :

1) Using AWK to select only reads that have a mapping quality of 60 like that :

samtools view -h sample1-bam1.bam | awk 'substr($0,1,1) == "@" || $5 == 60 {print}'| samtools view -bS - > sample1-MQ60-bam1.bam

bam2cfg.pl -h sample1-MQ60-bam1.bam sample1-MQ60-bam2.bam > sample1.cfg
breakdancer_max sample1.cfg > sample1-MQ60.txt  >> returns a file with found CNVs

2) In parallel on Breakdancer directly using -q flag

bam2cfg.pl -h sample1-1.bam sample1-2.bam > sample1-method2.cfg
breakdancer_max -q 60 sample1-method2.cfg > 1-MQ60.txt  >> returns an empty file
breakdancer_max -q 59 sample1-method2.cfg > 1-MQ59.txt  >> returns an empty file

Technically the two methods should return the same amount of CNV found .. why when I use -q flag on breakdancer this returns an empty file ??

Answer 1

[Chris Miller]

In the first code block, you're selecting reads with a mapping quality of exactly 60. In the second example, you're using reads that have a mapping quality of at least 60.

I'm not sure that completely explains your problem, but it's certainly somewhere to start. May I also ask why you're requiring a MQ of 60? That seems awfully high. 35 or 40 is a more typical cutoff for breakdancer, IIRC.

[madkitty]

I think by default it's 35 on breakdancer. But I had horrible long and multiple CNVs called whenever I use the default option on breakdancer, you can see it here : http://i46.tinypic.com/iwjdar.jpg

Even with mapping quality of 60 of the first method the result is not that different than the image you can see. If I use method 2) I have an empty file ..