[Bioc-devel] GenomicFiles: some random issues

Kasper Daniel Hansen kasperdanielhansen at gmail.com
Tue Sep 30 01:30:56 CEST 2014


Martin, I think Time explained it pretty well.  All I am saying is that it
would be convenient to have a class which represents that another class can
be treated as a GRanges.

In bsseq I have one such class: BSseqTstat.  I used to also construct the
BSseq class like this, but since it is essentially a SummarizedExperiment I
ended up changing the class definition, so it no longer has anything to do
with hasGRanges (which makes it a bit weird for a new reader; I used to
have hasGRanges to avoid code duplication between the two classes).

I have found this basic class very useful for setting up new classes,
indeed I have just used it in an internal package (which may be released
one day when it actually has useful stuff in it).  And I was pointing it
out in the context of GenomicFiles and SummarizedExperiment sharing
functionality but not sharing code through inheritance.

Michael: I am aware of additional columns in a GRanges object, but I
envision this in situations where I want to add extra slots to the basic
class.

Best,
Kasper



On Mon, Sep 29, 2014 at 5:43 PM, Michael Lawrence <lawrence.michael at gene.com
> wrote:

> I wrote it, because it might be clearer to see it as code:
>
> ###
> =========================================================================
> ### DelegatingGenomicRanges objects
> ###
> -------------------------------------------------------------------------
> ###
> ### Virtual class that delegates GenomicRanges data access to a
> ### GenomicRanges delegate.
> ###
>
> setClass("DelegatingGenomicRanges",
>          representation(delegate="GenomicRanges"),
>          contains=c("GenomicRanges", "VIRTUAL"))
>
> ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
> -
> ### Slot getters and setters.
> ###
>
> setMethod("seqnames", "DelegatingGenomicRanges",
>           function(x) seqnames(x at delegate))
>
> setMethod("ranges", "DelegatingGenomicRanges",
>           function(x, ...) ranges(x at delegate, ...))
>
> setMethod("strand", "DelegatingGenomicRanges", function(x)
> strand(x at delegate
> ))
> setMethod("seqinfo", "DelegatingGenomicRanges", function(x)
> seqinfo(x at delegate))
>
> ### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
> -
> ### Demo of using this.
> ###
>
> setClass("GenomicRangesWithStats",
>          representation(stats="matrix"),
>          contains="DelegatingGenomicRanges")
>
> setMethod("extraColumnSlotNames", "GenomicRangesWithStats", function(x)
> "stats")
>
>
> On Mon, Sep 29, 2014 at 2:29 PM, Michael Lawrence <michafla at gene.com>
> wrote:
>
> > This comes down to an inheritance vs. composition decision, but I hope
> > everyone is aware of the extraColumnSlots mechanism in GRanges that makes
> > it easy to add additional "column-shaped" slots to a subclass of GRanges.
> > It may be that hasGRanges could take a similar approach, except support
> > composition-based designs and thus complement (or maybe replace?) the
> > existing inheritance-based approach. Probably the composition approach
> > should at least provide something that inherits from the virtual
> > GenomicRanges class, because it would gain most of GRanges and only need
> to
> > implement the low-level data accessors, which would simply delegate to
> the
> > GRanges slot.
> >
> >
> >
> > On Mon, Sep 29, 2014 at 2:07 PM, Tim Triche, Jr. <tim.triche at gmail.com>
> > wrote:
> >
> >> Just to support Kasper's proposal, I see this as being broadly useful.
> >>
> >>
> >> It's really common to be dealing with an object (say, a
> >> SummarizedExperiment-derived doohickey, or a transcriptDb-looking
> whatzit,
> >> or "bumps", or "dmrcoutput") that has metadata which "looks like" a
> >> GRanges.  More often than not, a generic like
> >>
> >> ## standardGeneric for "granges" defined from package "GenomicRanges"
> >> granges <- function (x, use.mcols = FALSE, ...) { ... }
> >>
> >> gets defined in some derived class and the GenomicRanges generic
> >> specialized to the class.  Or sometimes because of S3-to-S4 impedance
> >> mismatch you instead get
> >>
> >> ## DMRcate results
> >> setClass('dmrcate.output')
> >> setMethod('granges', 'dmrcate.output', function(x, ...)
> extractRanges(x) )
> >>
> >> but the bottom line is that for all of these "things", if they contain a
> >> virtual class like "bsseq::hasGRanges", it enables semi-automatic
> >> extraction of their associated GRanges via the already-in-GenomicRanges
> >> granges() generic.  It seems handy to me, after writing some ungodly
> >> number
> >> of kludges to do similar things (see above).
> >>
> >> So if ::hasGranges could go into GenomicRanges or biobase or whatever as
> >> an
> >> abstract/virtual base class from which to inherit, with a mandatory
> method
> >> (is it possible to mandate a method for an S4 class?  been a while for
> >> me),
> >> then people who just want their intervals could reliably get them from
> >> output that hasGranges.
> >>
> >> ## what usually seems to happen
> >> foo <- emitWackyClassThatHasSomethingResemblingAGRangesAttachedToIt(x,
> y,
> >> z)
> >>
> >> ## the part that would be nice
> >> granges(foo) ## since the author made WackyClass contain "hasGranges"
> >>
> >> bsseq::hasGranges already defines a number of methods along these lines.
> >> More concretely, and stolen directly from bsseq/R/BSseqTstat_class.R:
> >>
> >> setClass("BSseqTstat", contains = "hasGRanges",
> >>          representation(stats = "matrix",
> >>                         parameters = "list")
> >>          )
> >> setValidity("BSseqTstat", function(object) {
> >>     msg <- NULL
> >>     if(length(object at gr) != nrow(object at stats))
> >>         msg <- c(msg, "length of 'gr' is different from the number of
> rows
> >> of 'stats'")
> >>     if(is.null(msg)) TRUE else msg
> >> })
> >>
> >> That's pretty handy and could stand to be in a very fundamental library
> >> IMHO.
> >>
> >>
> >>
> >> Statistics is the grammar of science.
> >> Karl Pearson <http://en.wikipedia.org/wiki/The_Grammar_of_Science>
> >>
> >> On Mon, Sep 29, 2014 at 12:58 PM, Martin Morgan <mtmorgan at fhcrc.org>
> >> wrote:
> >>
> >> > On 09/28/2014 07:49 PM, Kasper Daniel Hansen wrote:
> >> >
> >> >> 3. (well, this is GenomicRanges) It seems like GenomicFiles has a lot
> >> of
> >> >> similarities with SummarizedExperiment.  I have a lot of use cases
> for
> >> a
> >> >> simple class which has a single GRanges slot, and then support stuff
> >> like
> >> >> granges(), start(), end() etc.  Right now I have an implementation in
> >> >> bsseq:hasGRanges, but I think it could be widely useful and fits
> >> naturally
> >> >> with SummarizedExperiment.
> >> >>
> >> >
> >> > Hi Kasper -- not sure what you're trying to say here; sounds like a
> >> > GRanges (!) and bsseq::hasGRanges doesn't exist. Can you clarify
> please?
> >> >
> >> > Martin
> >> >
> >> > --
> >> > Computational Biology / Fred Hutchinson Cancer Research Center
> >> > 1100 Fairview Ave. N.
> >> > PO Box 19024 Seattle, WA 98109
> >> >
> >> > Location: Arnold Building M1 B861
> >> > Phone: (206) 667-2793
> >> >
> >> >
> >> > _______________________________________________
> >> > Bioc-devel at r-project.org mailing list
> >> > https://stat.ethz.ch/mailman/listinfo/bioc-devel
> >> >
> >>
> >>         [[alternative HTML version deleted]]
> >>
> >> _______________________________________________
> >> Bioc-devel at r-project.org mailing list
> >> https://stat.ethz.ch/mailman/listinfo/bioc-devel
> >>
> >
> >
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> Bioc-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/bioc-devel
>

	[[alternative HTML version deleted]]



More information about the Bioc-devel mailing list