Search
Question: How is the mean-variance plot calculated?
0
gravatar for neverstop
18 months ago by
neverstop0
neverstop0 wrote:

I'm trying to understand how the mean-variance plot is calculated with the meanVarPlot() of EDASeq package. I tried to reproduce it manually with the code below but the result is very different.

Given the expression matrix (row=genes, columns=samples), I used the logarithm to transform the counts before computing mean and variance. Since there are zero counts, I considered the  maximum between 0 and the log of the counts. After that I calculated the row means and the row variances (the mean and variance of each gene across all samples).

Here's the code I used:

rm(list=ls())
library(yeastRNASeq)
library(EDASeq)
data(geneLevelData)
head(geneLevelData)
Data=newSeqExpressionSet(as.matrix(geneLevelData))
meanVarPlot(Data,log=T)
geneLevelDataM=as.matrix(geneLevelData)

means=apply(pmax(log(geneLevelDataM),0),1,mean)
variances=apply(pmax(log(geneLevelDataM),0),1,var)
plot(means,variances)

means=apply(pmax(log2(geneLevelDataM),0),1,mean)
variances=apply(pmax(log2(geneLevelDataM),0),1,var)
plot(means,variances)

 

ADD COMMENTlink modified 18 months ago by Wolfgang Huber13k • written 18 months ago by neverstop0
2
gravatar for Wolfgang Huber
18 months ago by
EMBL European Molecular Biology Laboratory
Wolfgang Huber13k wrote:

Have a look at the source, e.g. https://github.com/Bioconductor-mirror/EDASeq/blob/master/R/methods-SeqExpressionSet.R (search for "meanVarPlot").

Other ways for accessing the source of a Bioconductor package are: download the .tar.gz from the package homepage and unpack; via the subversion server; via the R command line, you can list each function's body (it's a bit wordy for S4 methods).

It looks like the function does something else than what its manual page says, regarding the log argument.

 

ADD COMMENTlink modified 18 months ago • written 18 months ago by Wolfgang Huber13k

The function does indeed do something else than what the manual page says, i.e., it takes the log after computing mean and variance.

Thank you both for spotting this incongruence! I will fix it in devel.

Just to add to Wolfgang's answer re: source code, another way, within R, is to type:

getMethod("meanVarPlot", "SeqExpressionSet")

 

ADD REPLYlink written 18 months ago by davide risso490
Please log in to add an answer.

Help
Access

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