Question: Problem iterating through BAM file using pysam
1
gravatar for c9r1y
5.0 years ago by
c9r1y10
United States
c9r1y10 wrote:

Hello,

I am trying to use a for loop in python to loop through a bam file and take out reads that map to certain chromosomes, for further downstream analysis. I am having problems in the beginning, seemingly easiest part of my code.

samfile = pysam.Samfile("mybam.bam", "rb")

positions=["7:151970856", "7:151970856"]

for p in positions:
    start = p.split(":")
    chr=start[0]
    pos=start[1]
    chromreads=[]
    for read in samfile:
    if int(chr) == int(read.tid):
        chromreads.append(read)
    print(len(chromreads))

Output I'm getting:

2525527
0

Output I expect:

2525527
2525527

It seems like it's not even looping through the second position. What am I doing wrong?

Thank you.

pysam forloop bam sam python • 2.7k views
ADD COMMENTlink modified 2.5 years ago by Ram32k • written 5.0 years ago by c9r1y10
1
gravatar for c9r1y
5.0 years ago by
c9r1y10
United States
c9r1y10 wrote:

Answered my own question from StackOverflow

ADD COMMENTlink modified 2.5 years ago by Ram32k • written 5.0 years ago by c9r1y10
1
gravatar for Matt Shirley
5.0 years ago by
Matt Shirley9.5k
Cambridge, MA
Matt Shirley9.5k wrote:

You're iterating through the entire Samfile object, consuming it on the first pass of your inner for loop. When you arrive at the second iteration of your outer for loop, there is nothing left in the Samfile iterator and so it raises StopIteration which causes the inner loop to exit.

ADD COMMENTlink modified 2.5 years ago by Ram32k • written 5.0 years ago by Matt Shirley9.5k
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: 1602 users visited in the last hour
_