reads y-position in a bam viewer
1
0
Entering edit mode
6.7 years ago
sacha ★ 2.4k

Hi,

I m doing a BAM viewer. I wonder how to set read position in vertical axis ? I mean, given N reads , How can I display all of them to fit the space and allow smooth horizontal scrolling ?
So how to get y-coordinate given:

  • read width
  • read position
  • region viewport ( start, end )

I may need a pure function . Something like :

int ypos(region_start, region_end, read_size, read_begin ... )

But I m not sure.. See example in IGV bellow to understand

IGV

IGV bam position • 2.1k views
ADD COMMENT
2
Entering edit mode

A genuine question. What would be 1 thing about your proposed viewer that would make a compelling use case for it when we have IGV/IGB/tablet on one end and asciigenome at the other.

ADD REPLY
0
Entering edit mode

An application build with C++/Qt and seqan library under GPL license.

ADD REPLY
2
Entering edit mode

You could look at the source of other genome browsers (e.g. samtools tview, IGV, etc) to find out how they approached this problem. Also, search for pack rectangles, you will find plenty of resources - I think displaying reads on a genome browser is a particular case of rectangle packing.

ADD REPLY
1
Entering edit mode

It works ! Thanks guys for your help ! I just need some optimization now.

       QList<QList<seqan::BamAlignmentRecord> > rows;

        seqan::BamAlignmentRecord record;
        seqan::readRecord(record, bamFileIn);

        bool newRow = true;

        for (auto& row : rows)
        {
            seqan::BamAlignmentRecord rec = row.last();
            if (int(record.beginPos  > int(rec.beginPos + seqan::length(rec.seq))))
            {
                row.append(record);
                newRow = false;
                break;
            }
        }
        if (newRow)
        {
            QList<seqan::BamAlignmentRecord> row;
            row.append(record);
            rows.append(row);
        }

enter image description here

ADD REPLY
0
Entering edit mode

int(rec.beginPos + seqan::length(rec.seq)

don't use the read length, you'lll have to calculate the 'end' position using the cigar string.

ADD REPLY
0
Entering edit mode

Isn't that something you decide (based on what a user is doing in terms of clicking/selecting to zoom-in and out). You set the coordinates of the region user is looking at and that allows you to lookup how many reads there are in that windows and associated max depth? Not being a programmer perhaps I am missing something?

ADD REPLY
0
Entering edit mode

See the picture example bellow. I would like to know where to draw the red read. ( bottom or top) enter image description here

ADD REPLY
0
Entering edit mode

pack the rectangle has already answered that part. I assume you will also be providing a feature to "group reads according to read strand"? Then it may need to go at the bottom, if color represents read strand.

ADD REPLY
3
Entering edit mode
6.7 years ago

as said h.mon, pack the rectangle. One example of basic packing algorithm in my sources:

    List<List<SAMRecord>> rows = new ArrayList<>();
private void add(final SAMRecord rec)
    {
    // pileup
    for(final List<SAMRecord> row:this.rows)
        {
        final SAMRecord last=row.get(row.size()-1);
        if(right2pixel.apply(last)+ this.minDistanceBetweenPairs > left2pixel.apply(rec)) continue;
        row.add(rec);
        return;
        }

    final List<SAMRecord>  row=new ArrayList<>();
    row.add(rec);
    this.rows.add(row);
               }
ADD COMMENT
1
Entering edit mode

Ok thanks ! I assume I must use this recursive 'add' function in the following way :

>         for record in readBam(chromosom, start, end): 
>                add(record)
ADD REPLY

Login before adding your answer.

Traffic: 2798 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