[BioC] Coercion of MethyLumiSet to MethyLumiR changes raw methylation values?
Nick.Fishbane at hli.ubc.ca
Tue Mar 11 19:43:38 CET 2014
I am a data analyst currently experimenting with normalization methods for data from the Illumina 450K platform, and currently I'm examining what the different normalization packages (minfi, lumi, wateRmelon) offer, and how they differ from one another. In looking at methylumi/lumi, I ran into a question that I haven't been able to solve -- I've looked at the documentation, mailling list, and hacked as much of the internal functions as I am able -- how does lumiMethyR change the raw values from the input file? I noticed that when I extract values using "exprs", I get different raw methylation values compared to both the raw data and the data read by methylumiR, where my understanding was that lumiMethyR is just a wrapper for the former function.
Let me be explicit with an example... I'll use the raw data provided in the package:
> filename <- system.file('extdata/exampledata.samples.txt', package = 'methylumi')
The methylation values are based on the raw Illumina output. In this example, we look at the first probe reading from the first sample:
>raw_dat <- read.delim(filename)
For the very first probe ("AATK_E63_R"), we have an unmethylated signal ("CY3") of 1565.675, and a methylated signal ("CY5") 3054.818, so the beta calculated according to Illumina standards, M/(U+M+100), is displayed ("AVG_Beta") as 0.6471396.
Now when I use the function methylumiR to read the data, and convert M-values to beta-values (I use the function "m2beta" below), I get the same thing:
> m2beta <- function(x) 2^x/(1 + 2^x)
> mlsdat <- methylumiR(filename)
mlsdat is of the class "MethyLumiSet"
But when I use the function lumiMethyR, I get a different beta-value (e.g. 0.65445962 for the first probe)
> mlmdat <- lumiMethyR(filename)
mlmdat is of the class "MethyLumiM"
And I noticed that this change occurs in the second line of the "lumiMethyR" function, when the "MethyLumiSet" object is coerced to class "MethyLumiM".
>mlmdat <- as(mlsdat, "MethyLumiM")
So I guess my main question is: what is the coercion of class doing to the numerical calculation of methylation signals? Any other insight would be very much appreciated! And I apologize in advance if this question has already been addressed. For the record, I am using R version 3.0.2 on 64-bit Windows 7.
More information about the Bioconductor