gravatar for va90

2 hours ago by

Hi,
I am using pysam.pileup to get the base per information from a bam file.
I am using the following code:

import pysam
samfile = bamfh=pysam.AlignmentFile(bam_file,'rb')
def get_base_info(samfile,chrom,pos,map_qual):
for pileupcolumn in samfile.pileup( chrom, pos, pos+1, truncate= True):
    pileupcolumn.set_min_base_quality(0)
    if pileupcolumn.pos == pos:
        for pileupread in pileupcolumn.pileups:
            if pileupread.alignment.mapping_quality >= map_qual:
                if pileupread.alignment.is_reverse:
                    strand = "-"
                else:
                    strand = "+"
                read_len= pileupread.alignment.query_alignment_length
                qname= pileupread.alignment.query_name
                tag= pileupread.alignment.get_tag('RG')
                flag= pileupread.alignment.flag

                if pileupread.is_del and not pileupread.is_refskip:

                    print(flag,pileupread.alignment.mapping_quality,qname,"Deletion",pos)
                else:
                    qual= pileupread.alignment.query_qualities[pileupread.query_position]
                    read_base= pileupread.alignment.query_sequence[pileupread.query_position]
                    print(read_len,strand,flag,pileupread.alignment.mapping_quality,read_base,qual,qname,pos,tag,sep='t')
# for example for cordinate chr1 at position 100000 for reads with mapping quality higher than 20
get_base_info(samfile,chrom,100000,20)

by this code, I can get deletions, matches, and mismatches. Anyone knows how I can get insertions using pileup?

Thanks,
Vahid.

link

modified 26 minutes ago

written
2 hours ago
by

va9010



Source link