In insightsengineering/hermes: Preprocessing, analyzing, and reporting of RNA-seq data

This vignette explores what could be a good class inheritance design in the hermes package.

Objectives

We want to have two classes: - HermesData (short HD) which is a "normal" SummarizedExperiment (short SE) - RangedHermesData (short RHD) which should have row ranges in addition.

Both of the classes should fulfill same validation requirements, they need to have certain rowData and colData columns and specific assay ("counts").

Question

The question is now how these classes relate to each other, as well as if we want to inherit from RangedSummarizedExperiment (short RSE).

Idea 1: parallel inheritance

Here we inherit in parallel: - HD inherits from SE - RHD inherits from RSE

For illustration we work with the acronyms here so we don't intersect with the classes already in this package.

library(SummarizedExperiment)

HD <- setClass(
  "HD",
  contains = "SummarizedExperiment"
)

RHD <- setClass(
  "RHD",
  contains = "RangedSummarizedExperiment"
)

Validation

Let's see how we can define the validation methods here.

We can define a function:

.validate.HD <- function(object) {
  msg <- NULL
  if (!("counts" %in% assayNames(object))) {
    msg <- c(msg, "counts needs to be an assay")
  }
  if (!("A" %in% names(colData(object)))) {
    msg <- c(msg, "column A needs to be in colData")
  }
  if (!("B" %in% names(rowData(object)))) {
    msg <- c(msg, "column B needs to be in rowData")
  }
  msg
}

Then set the validity methods for both classes to this function:

setValidity2(
  "HD",
  method = .validate.HD
)

setValidity2(
  "RHD",
  method = .validate.HD
)

Then try this out:

se <- SummarizedExperiment(
  list(counts = matrix(1:4, 2, 2)),
  colData = data.frame(A = c(1, 2)),
  rowData = data.frame(B = c(3, 4))
)
hd <- HD(se)

This works, but:

RHD(se)

does not work. Seems we need to create the RSE differently.

rse <- as(se, "RangedSummarizedExperiment")
rhd <- RHD(rse)

OK so it works like that, nice.

Shared methods

How can we define the shared methods that work the same way, e.g. the counts accessor method, here?

One idea is to define a class union first.

setClassUnion(
  "AHD",
  c("HD", "RHD")
)

Then the define the method for that super class:

setMethod(
  "counts",
  signature = "AHD",
  definition = function(object) {
    assay(object, "counts")
  }
)

Try out:

counts(hd)
counts(rhd)

So that works well. Actually, since the validity method is the same, we can also do the same for the validation.

Coercion

We might also want to allow users to coerce between HD and RHD. We can walk the inheritance path up to the SE level and then use the coercion between SE and RSE before calling the constructor.

setAs(
  from = "RHD",
  to = "HD",
  def = function(from) {
    rse <- as(from, "RangedSummarizedExperiment")
    HD(as(rse, "SummarizedExperiment"))
  }
)

setAs(
  from = "HD",
  to = "RHD",
  def = function(from) {
    se <- as(from, "SummarizedExperiment")
    RHD(as(from, "RangedSummarizedExperiment"))
  }
)

Try that this works:

as(rhd, "HD")
as(hd, "RHD")

Row ranges

We will want to use the rowRanges() method on RHD (and HD where nothing is there).

rowRanges(rhd)
rowRanges(hd)

So this works. The big advantage here is that we get this for free via inheriting from RSE and we don't need to do anything ourselves.

Binding

rbind(rhd, rhd)
cbind(rhd, rhd)

So just works.

Constructor from (R)SE

Here we could have a single function and then decide based on RSE vs. SE input class.

library(assertthat)
MakeHD <- function(object) {
  assert_that(is_class(object, "SummarizedExperiment"))
  # check here things, annotate etc.
  if (is(object, "RangedSummarizedExperiment")) RHD(object) else HD(object)
}

Try out:

MakeHD(se)
MakeHD(rse)

Constructor from matrix

Also here we could have a single function.

HDfromMatrix <- function(counts, ...) {
  # se will either be SE or RSE, depending on ... contents.
  se <- SummarizedExperiment(
    list(counts = counts),
    ...
  )
  # Then MakeHD will do depending on that HD or RHD.
  MakeHD(se)
}

Try this out.

counts <- matrix(1:4, 2, 2)
hd <- HDfromMatrix(counts, 
                   colData = data.frame(A = c(1, 2)), 
                   rowData = data.frame(B = c(1, 2)))
hd

rowRanges <- GRanges(
  c("chr1", "chr2"),
  IRanges(c(124L, 134214L), width=100),
  strand=c("+", "-"),
  feature_id= c(1L, 2L),
  B = c(1, 2)
)
rhd <- HDfromMatrix(counts, rowRanges = rowRanges, colData = data.frame(A = c(1, 2)))
rhd