# filter secondary alignments using pysam

I am new to python and trying to learn. The below is an attempt to filter out secondary reads in a bam using pysam. I ran mark duplicates and got an error on several bam files due to secondary alignments. I added comments to each line to illustrate my thought process. Thank you :).

The specific mark duplicates error:

Value was put into PairInfoMap more than once


Here is the script:

#! /usr/bin/python  ## call python script
import sys ## import python system functions
import pysam  ## import module

bam = pysam.AlignmentFile(*.bam, "rb")  ## open bam and read
click.echo("Reading BAM file")  ## output message
for read in bam:  ## start loop and iterate over each bam
if read.is_secondary=true  ## not the primary alignment
if read.has_tag('XS')  ## look for Secondary Alignment Score tag
if read.has_tag('AS')  ## look for Alignment Score tag
if AS < XS  ## Alignment score less then
bam.close()  ## end processing

• Is there an issue with the script? If so what is it? Are you looking for feedback on the script? – Bioathlete Apr 13 at 2:34

Quick look at your script suggests that you are comparing the XS and AS from different reads. I would recommend removing the continue from the read.has_tag('XS') continue block.

I think that if you filter out is_secondary and is_supplementary you don't need to do the XS to AS comparison assuming whatever generated for BAM did so correctly.

for read in bam.fetch()
continue
if read.has_tag('XS') <- This block stores XS and then moves on to the next read

• I have not tried to executec it as this is my first python script i am looking for suggestions. Thank you :). – justaguy Apr 13 at 5:25
• I ran the script (modified it in the post) and there are many syntax errors... being new to python I am not sure what is best... what do you suggest? Thanks :). – justaguy Apr 15 at 13:12
• I would start will researching how to write to a file. read.write(read) is not correct syntax. You will need to open a file and write to that. Go error by error until you have a working script then determine if the out is correct. – Bioathlete Apr 15 at 16:03