Entering edit mode
4.5 years ago
va90
▴
50
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.
Hi @va90, I think you can do
This will reproduce mpileup format, from which you can infer all types of variants.
Another shortcut is to use
pileupread.indel
, it returns positive values incidcating insertions and insert lengh after this base and negative values indicating deletion legth after it.