Quick comments below: On Mon, Apr 6, 2009 at 11:38 AM, Angel <rubio.diaz.an...@gmail.com> wrote: > > Hi Group, > > I have been working with some expression arrays. We first did the > analysis in aroma.affymetrix, using RmaBackgroundCorrection, > quantileNormalization and RmaPlm. The arrays come from different > batches. This effect can be seen in the summarized expression. > > After the RmaBackgroundCorrection the plotDensity shows that there are > different group and the offset has not been completely removed. > plotDensity(csExprBC, type = “pm”)
The implement of RmaBackgroundCorrection utilizes bg.adjust() of the 'affy' package, which fits the norm+exp model and it always fits and corrects PM-only probes. > > We tried to use the function “backgroundCorrect” from the Limma > package using a maximum likelihood estimator for the parameters of the > normexp model. This is what we did: > > Ø aux <- getUnitIntensities(csExpr, stratifyBy="pm"); > Ø aux <- unlist(aux); > Ø dim(aux) <- c(41, length(aux)/41); # i have 41 samples > Ø ExprPM <- t(aux) ; It's better to do as follows: cdf <- getCdf(csR); cells <- getCellIndices(cdf, stratifyBy="pm", unlist=TRUE, useNames=FALSE); Y <- extractMatrix(csR, cells=cells); > > aux will be a matrix with whose size is genes (we used the ensembl cdf > for hgu133plus2) by samples, with only “pm” probes. Once I have the > variable “aux” I call the function from limma: > > Ø ExprBCmle <- backgroundCorrect(ExprPM, method="normexp", > normexp.method="mle",verbose=TRUE); > > > The density plot result is impressive: > > plotDensity(log2(ExprBCmle+50)) First thing to test is to see how much bg.adjust() of affy and backgroundCorrect() of limma differ. You might have to try different flavors of 'normexp.method', but "hopefully" you get something similar to what bg.adjust() (i.e. RmaBackgroundCorrection) gives. Make sure to compare without shifting 50 units! plotDensity(csExprBC, type=“pm”); plotDensity(log2(ExprBCmle)); If those are rather similar, then you know that you can use RmaBackgroundCorrection. Shifting probe signals with a constant is simple in downstream methods. Cheers, /Henrik > > > All the arrays seem to have almost identical densityplot (without > performing any quantile normalization). We summed up 50 as suggested > in (Gordon Smyth et al., Biostatistics January 2009. Microarray > background correction: maximum likelihood estimation for the normal- > exponential convolution). > > So my question is if there is anyway to do this from aroma.affymetrix, > I mean if this backgroundCorrection function is implemented in > aroma.affymetrix? If it is not, can we create a cel file with the > results of this background correction and continue afterwards? > This is the code we used: > >> sessionInfo() > > R version 2.8.1 (2008-12-22) > > i386-pc-mingw32 > locale: > LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States. > 1252;LC_MONETARY=English_United States. > 1252;LC_NUMERIC=C;LC_TIME=English_United States.1252 > > attached base packages: > > [1] stats graphics grDevices utils datasets methods > base > > other attached packages: > > [1] limma_2.16.4 aroma.affymetrix_1.0.0 > aroma.apd_0.1.3 > [4] R.huge_0.1.6 affxparser_1.14.0 > aroma.core_1.0.0 > [7] aroma.light_1.11.1 digest_0.3.1 > matrixStats_0.1.3 > [10] R.rsp_0.3.4 R.cache_0.1.7 > R.utils_1.1.3 > [13] R.oo_1.4.6 R.methodsS3_1.0.3 > > > # this is the code ive used, > > library(aroma.affymetrix) > verbose <- Arguments$getVerbose(-8) > timestampOn(verbose) > projectName <- "EndotelioExpr" > > chipType <- "HGU133Plus2_Hs_ENSG"; > cdf <- AffymetrixCdfFile$byChipType(chipType); > csExpr <- AffymetrixCelSet$byName(projectName,cdf=cdf); > > > # Background removal > ExprBC <- RmaBackgroundCorrection(csExpr); > ExprcsBC <- process(ExprBC, verbose = verbose); > plotDensity(ExprcsBC, type=”pm”); > > > # Background removal with limma > > aux <- getUnitIntensities(csExpr, stratifyBy="pm"); > aux <- unlist(aux); > dim(aux) <- c(41, length(aux)/41); > aux <- t(aux) ; > ExprBCmle <- backgroundCorrect(aux, method="normexp", > normexp.method="mle",verbose=TRUE); > plotDensity(log2(ExprBCmle+50)) > > > Thanks in advance for the answer. Best regards, > > Angel > > PD: I can post the density plots if anybody knows how to do it :-P > > > --~--~---------~--~----~------------~-------~--~----~ When reporting problems on aroma.affymetrix, make sure 1) to run the latest version of the package, 2) to report the output of sessionInfo() and traceback(), and 3) to post a complete code example. You received this message because you are subscribed to the Google Groups "aroma.affymetrix" group. To post to this group, send email to aroma-affymetrix@googlegroups.com To unsubscribe from this group, send email to aroma-affymetrix-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/aroma-affymetrix?hl=en -~----------~----~----~----~------~----~------~--~---
