Question: problem with getting a subset of alignments mapping to a GRange using readGAlignmentsFromBam
0
4.7 years ago by
rasi198310
United States
rasi198310 wrote:

I am unable to read a subset of alignments from a Bamfile as specified by a GRange. I am new to Bioconductor - am I missing something obvious?

779 Alignments in the Bamfile

 > readGAlignmentsFromBam(Bamfile) GAlignments object with 779 alignments and 0 metadata columns: seqnames strand cigar qwidth start end width njunc [1] 12 - 24M 24 25359522 25359545 24 0 [2] 12 - 30M 30 25359527 25359556 30 0 [3] 12 - 28M 28 25359755 25359782 28 0 [4] 12 - 27M 27 25362687 25362713 27 0 [5] 12 - 28M 28 25362734 25362761 28 0 ... ... ... ... ... ... ... ... ... [775] 12 - 22M5355N12M 34 25398308 25403696 5389 1 [776] 12 - 31M 31 25403736 25403766 31 0 [777] 12 - 25M 25 25403737 25403761 25 0 [778] 12 - 38M 38 25403745 25403782 38 0 [779] 12 - 24M 24 25403745 25403768 24 0 ------- seqinfo: 25 sequences from an unspecified genome  The my_grange object overlaps with only a fraction of the above alignments (92/779) > my_grange GRanges object with 1 range and 1 metadata column: seqnames ranges strand | tx_id | [1] 12 [25358180, 25362845] - | 47893,47894 ------- seqinfo: 1 sequence from an unspecified genome Using the which param after subseting with my_grange still returns all the 779 alignments > readGAlignmentsFromBam(Bamfile, which = ScanBamParam(which = my_grange)) GAlignments object with 779 alignments and 0 metadata columns: seqnames strand cigar qwidth start end width njunc [1] 12 - 24M 24 25359522 25359545 24 0 [2] 12 - 30M 30 25359527 25359556 30 0 [3] 12 - 28M 28 25359755 25359782 28 0 [4] 12 - 27M 27 25362687 25362713 27 0 [5] 12 - 28M 28 25362734 25362761 28 0 ... ... ... ... ... ... ... ... ... [775] 12 - 22M5355N12M 34 25398308 25403696 5389 1 [776] 12 - 31M 31 25403736 25403766 31 0 [777] 12 - 25M 25 25403737 25403761 25 0 [778] 12 - 38M 38 25403745 25403782 38 0 [779] 12 - 24M 24 25403745 25403768 24 0 ------- seqinfo: 25 sequences from an unspecified genome
ADD COMMENTlink
modified 4.7 years ago by James W. MacDonald52k • written 4.7 years ago by rasi198310
Answer: problem with getting a subset of alignments mapping to a GRange using readGAlign
2
4.7 years ago by
United States
James W. MacDonald52k wrote:

You want to pass the GRanges in using the 'param' argument, not the (non-existant) 'which' argument.

readGAlignmentsFromBam(Bamfile, param = ScanBamParam(which = my_grange))
ADD COMMENTlink written 4.7 years ago by James W. MacDonald52k

Thanks, I missed an R gotcha: you can pass nonexistent parameters and not get any errors!

ADD REPLYlink written 4.7 years ago by rasi198310
1

Well, that's because readGAlignmentsFromBam() has the ... argument. In BioC devel, readGAlignmentsFromBam() is deprecated in favor of readGAlignments(). And readGAlignments() doesn't have the ... argument anymore so you'll get an error if you try to pass a non existent argument to it. So it's not R but the way we design our functions and how we handle the ...

H.

ADD REPLYlink written 4.7 years ago by Hervé Pagès ♦♦ 14k

Thanks for the clarification.

ADD REPLYlink written 4.7 years ago by rasi198310
Please log in to add an answer.

Content
Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 291 users visited in the last hour