2
$\begingroup$

I have a BAM file (DNA, shallow whole genome sequencing at ~1X) where I want to filter reads (using sambamba) to keep only those which have a template length > 20 and mapping quality > 20, discarding unmapped reads.

The filter I have used is this:

sambamba view -h \
   -t 5 \
   -F "paired and template_length > 0 and mapping_quality > 20 and not (unmapped or mate_is_unmapped)" \
   -f bam 142-ready.bam \
   -o 142-filter.bam

However, inspection of the reads in IGV led me to think that perhaps reads on the minus strand are being somehow filtered out.

Pre filter:

Unfiltered view in IGV with mate pairs highlighted

Post filter:

Filtered view in IGV, showing only one mate pair

The only way I don't get such an aggressive filtering is by setting only not (unmapped or mate_is_unmapped). But that, obviously, defeats the purpose as it doesn't include what I need.

I can't, however, see what exactly is wrong with the filter. What is the correct way to filter for template length while not accidentally discarding reads on the minus strand?

$\endgroup$

1 Answer 1

1
$\begingroup$

The reason is that sambamba does not use the absolute value of the template length, so, if the read maps so the minus strand, the template length is actually negative. Setting a template length > 0 means that all reads mapping to the minus strand are removed.

In fact, you can't still filter for negative template lengths. See this GitHub PR.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.