Filter a SAM file by alignment score (AS)

[sjackman]

Hi,

I have a SAM file with the alignment score tag (AS) set appropriately (mapped by BWA-SW). I want to filter out all alignments with an alignment score (note: not mapping quality) less than a threshold. samtools doesn't seem to have this feature. Do you know of a software package that does?

Thanks,
Shaun

[swbarnes2]

You can probably just use grep or awk.

[nanto]

use:

samtools view -S file.sam | awk '{print $numberoffield}' > file

and just find which field is AS.

[Kulvait]

Solution

It is not so straightforward but for my bam file works well

samtools view -S file.bam | awk '{split($14, subfield, ":"); if(subfield[3]>120) print $0}', where 14 field match AS:i:NUM

[dpryan]

I would normally use pysam for something like this, just for the sake of simplicity.

[Kulvait]

Can you please provide full code you would use for filtering based on this paramter XS? This is sample row from my sam file to be used for filtering:

M03094:5:000000000-ALEAW:1:1102:21888:24571	16	chr9	139396897	60	27S213M27S	*	0	0	AAGCACAAAGGCGGCCGCCGCCACGTACATGAAGTGCAGCTGCGCTGGCGGGGGCGGCCCCCCGGGCTCACCTGCGGGCACGGGGGCCAGGGGCAGGGGCCCGGACATTAGGCAGCGGCTACGCAGCAGGCTGGTGGCCGGGGGGCGGCGGACTGGCTCCGCGTCCGGGCGCCTCCTCACCCACTCTCCTCCATCCCGCCCTCCAAAATAAGGTCATTCTCTACGCGATTAATCAGAATTGCAAACTATCGCTAAATTCTCTCCTGC	CCBBCFFFFFFCGGGGGGGGGGGGGHHHGHHHHHHHHHHHHHHGGGGGGGGGGG/<C--<A-:-----:0;F0/;-@F--..;[email protected];@-:-:-:BB0B009:FFF:FFFFFFFFFFFFFF9FG;FGGGGGG;GGGHHHHHGGGGGGGGGGGGGGGFGGGHHHHGGHGHHGHGGGHHGGGGGGGGGGHHGHHHHHHHHHHHHHHHHGGGGGGGHHHHHHHHHHHHHHHHHHHGGGGGGGGGGFFDCCDDDDDDC	MD:Z:45C12T2A3T31T10C109T48	NM:i:7	AS:i:232	XS:i:19

[dpryan]

import pysam
bam = pysam.AlignmentFile("foo.bam", "rb")
for read in bam.fetch():
    if read.has_tag('XS'):
        XS = read.get_tag('XS')
        # do something
    else:
        # do something else
bam.close()