On Wed, Nov 27, 2013 at 4:28 AM, Nicolas Delhomme <nicolas.delhomme@umu.se>wrote: > Hej Hervé! > --------------------------------------------------------------- > Nicolas Delhomme > > Nathaniel Street Lab > Department of Plant Physiology > Umeå Plant Science Center > > Tel: +46 90 786 7989 > Email: nicolas.delhomme@plantphys.umu.se > SLU - Umeå universitet > Umeå S-901 87 Sweden > --------------------------------------------------------------- > > On 22 Nov 2013, at 22:03, Hervé Pagès <hpages@fhcrc.org> wrote: > > > Hi Nico, > > > > On 11/22/2013 06:36 AM, Nicolas Delhomme wrote: > >> Hej Michael, Hervé! > >> > >> I have to admit that I’m not exactly following the details of your > discussion here. What would appear intuitive for me for a findOverlaps > function on GRangesList, GRangesList is that it returns a list of > findOverlaps on GRanges,GRanges. I can see difficulties in determining the > concordance between the GRanges of the two GRangesList though. And I don’t > have the overview to decide whether such a behaviour would be consistent, > but that’s naively what I would expect. > > > > Typically when 'query' and 'subject' are GRangesList, one is > > representing alignments (e.g. single-end with gaps or paired-end, > > thus more than 1 range is needed to represent an alignment), and > > the other one is representing a set of compound genomic features > > e.g. exons grouped by transcript or by gene. > > findOverlaps() then just returns yes/no answers to the questions > > "do query[[i]] and subject[[j]] overlap?" for all the possible (i, j) > > combinations. It doesn't return any detailed information about which > > ranges in GRanges objects query[[i]] and subject[[j]] have overlaps. > > If a hit is reported between query[[i]] and subject[[j]], it just > > means that at least one range in the former overlaps with one range > > in the latter. For many analyses, this is enough information i.e. > > knowing that there is a hit between a read and a transcript is enough. > > > > Sure, completely make senses. > > > For analyses that need to know a little bit more about the > > characteristics of those hits, some post-processing of the Hits object > > is generally required. That's what tools like summarizeOverlaps or > > findSpliceOverlaps do. In your case, if you want a "list of > > findOverlaps on GRanges,GRanges", it's also something that you could > > obtained by post-processing the hit object: > > > > ov <- findOverlaps(query, subject, ...) > > mapply(findOverlaps, query[queryHits(ov)], subject[subjectHits(ov)]) > > > > Exactly what I ended up doing and as my GRangesList are not large ~50 > elements, it went fast. > mclapply helped made it even faster obviously :-) > > > > It probably won't be very efficient but we could address this by > > implementing a fast "parallel findOverlaps" (e.g. pfindOverlaps) > > if there is interest in something like that. > > > > I have no such needs yet, but thanks. > > > But back to your original request of supporting type="equal" on > > GRangesList objects. When you use the 'type' argument, that only > > changes the criteria that findOverlaps uses for deciding whether > > to report a hit or not. That doesn't change what is reported. > > What is reported is still a Hits object answering the yes/no > > answers. When type="equal" will be available, a hit between > > query[[i]] and subject[[j]] will be reported only if the ranges > > are the same in the 2 objects. In that case there doesn't seem to > > be a lot of value in trying to generate a "list of findOverlaps on > > GRanges,GRanges" because you already know what you need to know. > > > > And that’s the behaviour I should indeed expect. I was biased by my own > needs ;-) > > > So what would be more useful: (1) a fast pfindOverlaps, (2) supporting > > type="equal" on GRangesList, (3) both? Knowing a little bit more > > about what you are trying to do would help us decide. > > > > I think (2) would be good simply because it is stated as possible in the > man page. Or the man page would need correction. > It's actually not stated in the man page, in fact the findOverlaps functionality on GRangesList is not well documented; there is no usage for the method, for example. Another useful thing along these lines would be == for GRanges(List). I and others have faced the use case of performing an overlaps operation, aligning the query and subject by indexing with query/subjectHits, and then wanting to know which were equal and which weren't. For example, when looking for compatibility between reads and transcripts. Metadata would be ignored, of course. > > > FWIW let me mention another way of post-processing a Hits object. > > encodeOverlaps() can be used on the Hits object to produce "overlap > > encodings". These are small strings that describe how the individual > > ranges in the 2 GRanges objects involved in a hit are positioned > > with respect to each other. The idea is a little bit the same as > > with the CIGAR encodings that provide details about how reads align > > to the reference. See the OverlapEncoding vignette in the GenomicRanges > > package for more info (it's a low level tool and it can be a little > > bit tricky to use, also the vignette is probably too technical and > > would need to be improved). > > That sounds very interesting; I’ll take a look at it. > > Thanks for this very detailed very useful answer! > > Cheers, > > Nico > > > > > Cheers, > > H. > > > >> > >> Anyway, just my 2 cts. Thanks for taking a look. > >> > >> Cheers, > >> > >> Nico > >> > >> --------------------------------------------------------------- > >> Nicolas Delhomme > >> > >> Nathaniel Street Lab > >> Department of Plant Physiology > >> Umeå Plant Science Center > >> > >> Tel: +46 90 786 7989 > >> Email: nicolas.delhomme@plantphys.umu.se > >> SLU - Umeå universitet > >> Umeå S-901 87 Sweden > >> --------------------------------------------------------------- > >> > >> On 21 Nov 2013, at 19:43, Michael Lawrence <lawrence.michael@gene.com> > wrote: > >> > >>> So I've checked into devel a match,GRangesList,GRangesList. This > allows findMatches() to return what you want. There is a question though > before this is approved: does it make sense for match() to act like > findOverlaps and consider each GRanges atomically (one returned index per > GRanges) or should match behave as it does other Lists and return an > IntegerList, with a value per range, grouped by the top-level elements. If > we decide on the latter, then the method I wrote needs to be removed and > the implementation moved to the "equals" mode in findOverlaps. Either way, > findOverlaps(type="equals") should be made to work. > >>> > >>> Michael > >>> > >>> > >>> On Thu, Nov 21, 2013 at 8:13 AM, Nicolas Delhomme < > nicolas.delhomme@umu.se> wrote: > >>> Thanks! > >>> --------------------------------------------------------------- > >>> Nicolas Delhomme > >>> > >>> Nathaniel Street Lab > >>> Department of Plant Physiology > >>> Umeå Plant Science Center > >>> > >>> Tel: +46 90 786 7989 > >>> Email: nicolas.delhomme@plantphys.umu.se > >>> SLU - Umeå universitet > >>> Umeå S-901 87 Sweden > >>> --------------------------------------------------------------- > >>> > >>> On 21 Nov 2013, at 17:06, Michael Lawrence <lawrence.michael@gene.com> > wrote: > >>> > >>>> I will work on this today. > >>>> > >>>> Michael > >>>> > >>>> > >>>> On Thu, Nov 21, 2013 at 4:43 AM, Nicolas Delhomme < > nicolas.delhomme@umu.se> wrote: > >>>> Hej Bioc! > >>>> > >>>> When I try to find “equal” ranges from two GRangesList object, I get > the following error: > >>>> > >>>>> findOverlaps(query=grng.def,subject=grng.mod,type="equal") > >>>> Error in match.arg(type) : > >>>> 'arg' should be one of “any”, “start”, “end”, “within” > >>>> > >>>> Isn’t type=“equal” supported for the GRangesList, GRangesList > signature? > >>>> > >>>> Cheers, > >>>> > >>>> Nico > >>>> > >>>> sessionInfo() > >>>> R version 3.0.2 (2013-09-25) > >>>> Platform: x86_64-apple-darwin13.0.0 (64-bit) > >>>> > >>>> locale: > >>>> [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8 > >>>> > >>>> attached base packages: > >>>> [1] parallel stats graphics grDevices utils datasets > methods base > >>>> > >>>> other attached packages: > >>>> [1] easyRNASeq_1.8.2 ShortRead_1.20.0 Rsamtools_1.14.1 > GenomicRanges_1.14.3 DESeq_1.14.0 lattice_0.20-24 > locfit_1.5-9.1 > >>>> [8] Biostrings_2.30.1 XVector_0.2.0 IRanges_1.20.5 > edgeR_3.4.0 limma_3.18.3 biomaRt_2.18.0 > Biobase_2.22.0 > >>>> [15] genomeIntervals_1.18.0 BiocGenerics_0.8.0 intervals_0.14.0 > >>>> > >>>> loaded via a namespace (and not attached): > >>>> [1] annotate_1.40.0 AnnotationDbi_1.24.0 bitops_1.0-6 > DBI_0.2-7 genefilter_1.44.0 geneplotter_1.40.0 grid_3.0.2 > hwriter_1.3 > >>>> [9] latticeExtra_0.6-26 LSD_2.5 RColorBrewer_1.0-5 > RCurl_1.95-4.1 RSQLite_0.11.4 splines_3.0.2 stats4_3.0.2 > survival_2.37-4 > >>>> [17] tools_3.0.2 XML_3.98-1.1 xtable_1.7-1 > zlibbioc_1.8.0 > >>>> > >>>> > >>>> --------------------------------------------------------------- > >>>> Nicolas Delhomme > >>>> > >>>> Nathaniel Street Lab > >>>> Department of Plant Physiology > >>>> Umeå Plant Science Center > >>>> > >>>> Tel: +46 90 786 7989 > >>>> Email: nicolas.delhomme@plantphys.umu.se > >>>> SLU - Umeå universitet > >>>> Umeå S-901 87 Sweden > >>>> > >>>> _______________________________________________ > >>>> Bioconductor mailing list > >>>> Bioconductor@r-project.org > >>>> https://stat.ethz.ch/mailman/listinfo/bioconductor > >>>> Search the archives: > http://news.gmane.org/gmane.science.biology.informatics.conductor > >>>> > >>> > >>> > >> > >> _______________________________________________ > >> Bioconductor mailing list > >> Bioconductor@r-project.org > >> https://stat.ethz.ch/mailman/listinfo/bioconductor > >> Search the archives: > http://news.gmane.org/gmane.science.biology.informatics.conductor > >> > > > > -- > > Hervé Pagès > > > > Program in Computational Biology > > Division of Public Health Sciences > > Fred Hutchinson Cancer Research Center > > 1100 Fairview Ave. N, M1-B514 > > P.O. Box 19024 > > Seattle, WA 98109-1024 > > > > E-mail: hpages@fhcrc.org > > Phone: (206) 667-5791 > > Fax: (206) 667-1319 > > [[alternative HTML version deleted]]