Search code examples
rlimma

How to get the underling data of a MA plot?

I would like to use the plotMA function of limma.

The example of the documentation works fine:

A <- runif(1000,4,16)
y <- A + matrix(rnorm(1000*3,sd=0.2),1000,3)
status <- rep(c(0,-1,1),c(950,40,10))
y[,1] <- y[,1] + status
plotMA(y, array=1, status=status, values=c(-1,1), hl.col=c("blue","red"))

Now I would like to access the underlying data that is used for the plot as I would like to use the data in a different context, not just the plot. I currently don't see a way to access the data; of course I could implement the method myself and only use the data, but that feels wrong.

Is there a way to access the underlying data used for the MA plot?

Solution

  • Looking at the code of plotMA we see that several variables are created and used for plotting. These variables are not returned however.

    You could now copy and paste the function to write your own function, which plots and returns the data. This is however, error-prone,if there is a new version of the function you may rely on old code.

    So what you can do instead is to use trace to insert arbitrary code into plotMA notably some code which stores the data in your global environment. I illustrate the idea with a toy example:

    f <- function(x) {
   y <- x + rnorm(length(x))
   plot(x, y)
   invisible()
}

    If we would like to use y in this function we could do something like this

    trace(f, exit = quote(my_y <<- y))
# [1] "f"
ls()
# [1] "f"
f(1:10)
# Tracing f(1:10) on exit 
ls()
# [1] "f"    "my_y"

    And now we can access my_y.

    What you should do:

    1. Look at the code of plotMA
    2. Identify which part of the data you need (e.g. x, y and sel)
    3. Use trace(plotMA, exit = quote({my_data <<- list(x, y, sel)}), where = asNamespace("limma"))
    4. Run plotMA
    5. Access the data via my_data

    Note. Check out ?trace to fully understand the possibilities of it. In particular, if you want to inject your code not at the end (exit) but at another psoition (maybe because intermediate variables are overwritten and you need the first results) for which you would need to use the at parameter of trace

    Update

    Maybe the easiest is to get a full dump of all local variables defined in the function:

    trace("plotMA", exit = quote(var_dump <<- mget(ls())), where = asNamespace("limma"))