htsjdk insert informations from "concordantly mapped" read pairs
1
0
Entering edit mode
7.0 years ago

Hi, I hope i did not post this in the wrong subforum or category.

I have a problem with using picard-htsjdk in a java program. I have alignments with paired-end sequencing data created in bowtie2, and i want to infer information about inserts between "concordantly" mapped read pairs.

I am using the following code to obtain the information from the SAMRecords:

        File bamFile = new File(this.alignmentPath);
        SamReaderFactory samReaderFactory = SamReaderFactory.makeDefault()
                .validationStringency(ValidationStringency.LENIENT);
        SamReader samReader = samReaderFactory.open(bamFile);


        for (final SAMRecord rec : samReader)
            if (rec.getReadPairedFlag())
                if (rec.getProperPairFlag())
                    if (rec.getFirstOfPairFlag()) {
                        System.out.println("First Start: " + rec.getAlignmentStart());
                        System.out.println("First End: " + rec.getAlignmentEnd());
                        System.out.println("Insert size: " + rec.getInferredInsertSize());
                        System.out.println("Mate End: " + rec.getMateAlignmentStart());
                        System.out.println();
                    }
        }

What I'm getting is something like this:

First Start: 64973835 First End: 64973887 Insert size: 75 Mate End: 64973829

First Start: 64976212 First End: 64976262 Insert size: -80 Mate End: 64976182

First Start: 64978477 First End: 64978522 Insert size: 129 Mate End: 64978572

First Start: 64 978 546 First End: 64 978 615 Insert size: 75 Mate End: 64 978 581

First Start: 64978610 First End: 64978662 Insert size: -129 Mate End: 64978533

First Start: 64978661 First End: 64978715 Insert size: -101 Mate End: 64978614

First Start: 64984272 First End: 64984313 Insert size: 79 Mate End: 64984317

First Start: 64 990 218 First End: 64 990 277 Insert size: -108 Mate End: 64990169

The Mate End seem not to make any sense at all, given that the information about the insert size and position of the "First in Pair" read is correct.

How can this be? Am I using the library in a wrong way, or is there something wrong with the alignment itself?

Thank you very much in advance for your help and opinions. Cheers! Stefan

alignment • 1.6k views
ADD COMMENT
0
Entering edit mode

last time I've checked, there was a problem with bowtie/tophap flags: Tophat : Sam-Flag 115 = Properly-Paired + Read.Reverse + Mate.Reverse ?

ADD REPLY
0
Entering edit mode

Oh, this is just a Mistake in the Println(). There does not seem to exist a getMateAlignmentEnd() Function.

ADD REPLY
0
Entering edit mode
7.0 years ago

Watch out you have "Mate End: " but the code says rec.getMateAlignment*Start*() shouldn't it be rec.getMateAlignmentEnd()? Maybe try also to look at some read pairs in a genome browser to see where these numbers come from.

ADD COMMENT

Login before adding your answer.

Traffic: 1839 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6