I have a BAM file:
@SQ SN:chr1 LN:248956422
@SQ SN:chrx LN:248956423
ST-E00110:348:HGVKKALXX:1:1201:5822:48670 323 chr1 9999 0 67H66M16H chrx 1000 0 GATAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAACCCTAAC JJJJJJJJJJJJJJJJAJJJJJJJJJJJJFJJJJJJFJFJJJJJJFJJJJJJJJJJJA77FJFJJJ NM:i:0 MD:Z:66 AS:i:66 XS:i:65 SA:Z:chr5,18606834,-,73S76M,34,0; RG:Z:g1
There is a read aligned to chr1
, and it's mate aligned to chrx
.
I have a BED file:
chr1 0 100000 TestOnly
I would like to filter out everything that falls outside my BED region, that includes cross-alignments. In my example, although my read aligned to chr1
but it's mate is not. I don't want this read.
When I do:
samtools view -L test.bed test.bam
the command gives me the read because it doesn't check cross-alignments.
My solution:
samtools view -L test.bed test.bam | grep -v chrx
but this is very slow and clumsy. In my production pipeline I would have to do something like:
samtools view -L test.bed test.bam | grep -v chrx | grep -v ... | grep -v ... | grep -v ... | grep -v ...
Q: Is there a better solution?