Question: pysam doesn't calculate correct nucleotide quantities from .bam file
1
gravatar for vit.filippo
26 days ago by
vit.filippo50
vit.filippo50 wrote:

I am very new to pysam, but I am trying to generate a table in which are listed all the nucleotides by genomic position divided by strand with python. I am working on amplicon data, so it would be possible that pysam would mark them as duplicates (as default). By now I simply re-wrote codes from a previous post (https://www.biostars.org/p/215192/#217314), but counts are not correct. In particular, I observed that only nucleotides observed in Read1 are reported (and are not fully correct), but counts from Read2 are almost absent, even if there are some positions with high counts (but presumably wrong). Anyone experienced such a strange behaviour? Below is reported the code for simplicity: import pysam

samfile = pysam.Samfile("reads.sorted.bam", "rb")
mystring = ''
mylist = []
for pileupcolumn in samfile.pileup("gi|251831106|ref|NC_012920.1|"):        
    mystring = ''
    for pileupread in pileupcolumn.pileups:
        if not pileupread.is_del and not pileupread.is_refskip:                
            mystring = mystring + pileupread.alignment.query_sequence[pileupread.query_position]
            mystring = "".join(sorted(list(mystring)))  # Sort alphabetically
            mylist.append((pileupcolumn.pos+1,mystring))
samfile.close()
python pysam count • 88 views
ADD COMMENTlink modified 26 days ago by dariober11k • written 26 days ago by vit.filippo50
0
gravatar for dariober
26 days ago by
dariober11k
WCIP | Glasgow | UK
dariober11k wrote:

Some random thoughts looking at the documentation for pysam/pileup:

  • max_depth: default to 8000, maybe you need to increase this limit

  • ignore_overlaps: the default is True and it may explain why only read 1 is counted

  • flag_filter: Maybe you want to disable this filter and count all reads

ADD COMMENTlink written 26 days ago by dariober11k

It worked like a charm. Thanks a lot, I didn't this section. Bests

ADD REPLYlink written 26 days ago by vit.filippo50
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1379 users visited in the last hour