Could not reproduce the results (MA-plot) in the paper when analyzing two-color Microarray data with Limma.

Question

[liux.bio]

Hi, biostars,

I am using Limma package to analyze a chip-on-chip data set generated by Mouse PromoterChip 5A.0 arrays. The arrays are two-color Microarray. The protocol the author used are as follows:

Extracted raw data were filtered and normalized using "normexp" method from the Limma package developed within the Bioconductor project in the R statistical programming environment (Smyth 2004). The two channels were balanced by lowess normalization using 0.3 as the span parameter with reduced weights in control and poor quality spots, followed by scaling across chips. The spots lower than 2-fold the local background in the both channels were removed from the study. The ratio of expression for each element on the array was calculated in terms of M (log2(Cy5/Cy3)) and A ((log2(Cy5) + log2(Cy3))/2)).

I followed the protocol and here is my code:

# load raw data 
library("limma")
targets <- readTargets("./rawData/E-MEXP-1714/raw/Targets.txt")
# filter 
filterFun <- function(x) { 
   flag <- as.numeric(x[,"GenePix:Flags"]) > -50.5
   okGreen <- as.numeric(x[,"GenePix:F532 Median"])/as.numeric(x[,"GenePix:B532 Median"]) > 2
   okRed <- as.numeric(x[,"GenePix:F635 Median"])/as.numeric(x[,"GenePix:B635 Median"]) > 2

   #as.numeric(flag && (okGreen || okRed)) 
   return(as.numeric(flag & (okGreen | okRed)))
 }

RG<-read.maimages(paste("./rawData/E-MEXP-1714/raw/",targets$FileName, sep=""), 
              columns = list(Gf="GenePix:F532 Median", Gb="GenePix:B532 Median", 
                           Rf="GenePix:F635 Median", Rb="GenePix:B635 Median"),
             annotation = c("Reporter identifier", "GenePix:Flags"),
              wt.fun = filterFun)


RG <- backgroundCorrect(RG, method="normexp")
RG <- normalizeWithinArrays(RG, method="loess")
MA <- normalizeBetweenArrays(RG, method = "Aquantile")

plotMA(MA, xlab="A(log2(ChIP + input)/2)", ylab="M(log2ChIP/input)", main="")

However, the MA-plot I plot is different from the one in the paper. The first figure is from the paper, the second one is the one I plot.

I am a newbie for Microarray data analysis. Could you please tell me what I missed?

Answer 1

[Gordon Smyth]

It can be hard to reproduce published analyses exactly, unless the authors provide exact code, but I can see several discrepancies between the authors' stated protocol and your code.

The most important discrepancy is that you haven't filtered any spots. The authors removed spots with low intensities whereas you have kept every spot in your plot.

The wt.fun argument is for defining weights rather than for filtering. The authors said that they defined reduced weights for poor quality spots, so you need to figure out what the reduced weight was (it apparently was not zero) and how they defined poor quality weights (it might have been in terms of Flag but maybe not).

The filtering needs to be done separately to computing the weights. It isn't clear at what stage the low intensity spots were removed. Probably after the background correction and normalization steps.

The authors do not mention having done between array normalization.

You could try adding zero.weights=FALSE to your plotMA function call to remove the spots with zero weight from the plot. But it still won't be quite the same as the published analysis unless you can attend to the above points.

[liux.bio]

Many thanks for your clear answer. I will have a try.