Question: bam stats F1F2/F1R2/R1F2/R1R2/F2F1/F2R1/R2F1/R2R1
2
gravatar for 14134125465346445
2.8 years ago by
United Kingdom
141341254653464453.5k wrote:

Hi,

I would like to know the statistics of the different read pair orientations for a bam file aligned with bwa mem.

I tried bamtools stats, but it does not split the counts into all possible orientations: F1F2/F1R2/R1F2/R1R2/F2F1/F2R1/R2F1/R2R1

Is there a piece of software that will give me the stats of each orientation possible?

Thx

bam • 2.0k views
ADD COMMENTlink modified 2.8 years ago by John12k • written 2.8 years ago by 141341254653464453.5k
1
gravatar for Pierre Lindenbaum
2.8 years ago by
France/Nantes/Institut du Thorax - INSERM UMR1087
Pierre Lindenbaum124k wrote:

using bioalcidae: https://github.com/lindenb/jvarkit/wiki/BioAlcidae

var count={};

while(iter.hasNext())
    {
    var rec = iter.next();
    if(!rec.getReadPairedFlag()) continue;
    if(rec.getReadUnmappedFlag()) continue;
    if(rec.getMateUnmappedFlag()) continue;
    /* I'm not sure you want this */
    if(rec.getAlignmentStart() >= rec.getMateAlignmentStart()) continue;

    if(rec.getReadNegativeStrandFlag())
        {
        key="R";
        }
    else
        {
        key="F";
        }

    key+= (rec.getFirstOfPairFlag()?"1":"2");

    if(rec.getMateNegativeStrandFlag())
        {
        key+="R";
        }
    else
        {
        key+="F";
        }

    key+=(rec.getFirstOfPairFlag()?"2":"1");

    if(!(key in count))
        {
        count[key]=1;
        } 
    else
        {
        count[key]++;
        }
    }

for(var key in count)
    {
    out.println(key+" "+count[key]);
    }

usage:

 java -jar src/jvarkit/dist/bioalcidae.jar -f file.js in.bam

or

 samtools view -F (something) -bu in.bam | java -jar src/jvarkit/dist/bioalcidae.jar -f file.js -F BAM
ADD COMMENTlink modified 2.8 years ago • written 2.8 years ago by Pierre Lindenbaum124k
1
gravatar for John
2.8 years ago by
John12k
Germany
John12k wrote:

I hate competing with Pierre because I always lose :P but SeQC can do this so I ought to plug it.

Download all the code from https://github.com/JohnLonginotto/SeQC You want to use the TYPE module, which will give a number for all the possible orientations a single or paired end read can align, as a number from 1 to 20. To decode this number, use this chart: http://i.imgur.com/DCMrtwe.jpg

If the image quality isn't good enough or you want help running SeQC, i'll be sitting in ac.gt/chat for the next 4 hours to help

ADD COMMENTlink written 2.8 years ago by John12k
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: 1403 users visited in the last hour