Question: Samlocusiterator Skips Every Second Locus
3
gravatar for Johan
6.2 years ago by
Johan840
Sweden
Johan840 wrote:

I'm writing some Scala code which uses the SamLocusIterator after a excellent tip here: A: Java fasta reference iterator

However I'm getting some pretty weird output, which makes it seem like the iterator is skipping every other locus. Here's an example:

pos: chr1:99984
pos: chr1:99986
pos: chr1:99988
pos: chr1:99990
pos: chr1:99992
pos: chr1:99994
pos: chr1:99996
pos: chr1:99998
pos: chr1:100000

I can replicate this both in Scala and in Java with the following example programs:

The java version

import java.io.File;
import java.util.Iterator;

import net.sf.picard.util.SamLocusIterator;
import net.sf.picard.util.SamLocusIterator.LocusInfo;
import net.sf.samtools.SAMFileReader;


public class JavaMinimalIteratorTest {
    /**
     * @param args
     */
    public static void main(String[] args) {
        File testFile = new File("exampleBAM.bam");
        SAMFileReader fileReader = new SAMFileReader(testFile);

        SamLocusIterator samLocusIterator = new SamLocusIterator(fileReader);
        Iterator<LocusInfo> iterator = samLocusIterator.iterator();

        while (iterator.hasNext()) {
            LocusInfo locus = iterator.next();
            System.out.println("pos: " + locus.getSequenceName() + ":" + locus.getPosition());
        }
    }
}

The Scala version

import net.sf.samtools.SAMFileReader
import net.sf.picard.util.SamLocusIterator
import java.io.File

object MinimalIteratorExample extends App {

    val testFile = new File("exampleBAM.bam")
    val fileReader = new SAMFileReader(testFile)

    val samLocusIterator: SamLocusIterator = new SamLocusIterator(fileReader)
    val iterator = samLocusIterator.iterator()

    while (iterator.hasNext()) {
        val locus = iterator.next()
        println("pos: " + locus.getSequenceName() + ":" + locus.getPosition())
    }
}

Am I'm using the iterator in the wrong way, or is this some bug in the Picard code? I've been looking debugging and looking into the Picard classes, but this far I've not been able to make sense of it. Any help figuring this out would be much appreciated.

java picard • 1.9k views
ADD COMMENTlink modified 5.1 years ago by agolicz30 • written 6.2 years ago by Johan840
3

I'm afraid 'hasNext()' calls 'next()' ...

ADD REPLYlink written 6.2 years ago by Pierre Lindenbaum119k

Well, I see. That's not great, now, is it? I'll post it as an issue to the picard/samtools-dev list.

ADD REPLYlink written 6.2 years ago by Johan840

Thank you. I just had a glance at github: the only place I saw this API used was a test in GATK (a simple loop to get the time for scanning a bam.... )

ADD REPLYlink written 6.2 years ago by Pierre Lindenbaum119k

did you post your question on the samtools mailing list ? I cannot see the post. http://sourceforge.net/mailarchive/forum.php?forum_name=samtools-devel

ADD REPLYlink written 6.2 years ago by Pierre Lindenbaum119k

I did, but only realized this morning that I had to be a member on the list to post. I posted it now though.

ADD REPLYlink written 6.2 years ago by Johan840

I think it is not the case. I just tried the similar code above replacing the "hasNext()" condition to true. But got the same result.

ADD REPLYlink written 6.2 years ago by Sangwoo Kim340

I've got the same problem here. Thanks I've never noticed this problem. you'd better ask the picard/samtools-dev mailing list.

ADD REPLYlink written 6.2 years ago by Pierre Lindenbaum119k

In case you did not know already: You can import scala.collection.JavaConversions.asScalaIterator and you the SamLocusIterator as every iterator. e.g.: for( entry <- iterator) { do something meaningful with entry } Saves some LOC and looks much cleaner.

ADD REPLYlink modified 4.7 years ago • written 4.7 years ago by Fabian Bull1.3k
2
gravatar for Sangwoo Kim
6.2 years ago by
Sangwoo Kim340
UC San Diego
Sangwoo Kim340 wrote:

I just found that this occurs in regions where no reads have been aligned. Because these regions are not of our interest, we should skip them by setting setEmitUncoveredLoci()

Try this.

SamLocusIterator samLocusIterator = new SamLocusIterator(fileReader);
samLocusIterator.setEmitUncoveredLoci(false);
Iterator<LocusInfo> iterator = samLocusIterator.iterator();
ADD COMMENTlink written 6.2 years ago by Sangwoo Kim340
1

You're right, this only happens in uncovered regions. :) Still somewhat unexpected, though.

ADD REPLYlink written 6.1 years ago by Johan840

Nice to know ! Thank you.

ADD REPLYlink modified 6.1 years ago • written 6.1 years ago by Pierre Lindenbaum119k
0
gravatar for agolicz
5.1 years ago by
agolicz30
Australia
agolicz30 wrote:

This bug is fixed now, as of picard tools v1.88. In case somebody also needs non-covered regions.

ADD COMMENTlink modified 21 months ago • written 5.1 years ago by agolicz30
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: 645 users visited in the last hour