[BioC] limma Normalization question

[Cecilia McGregor]

Dr Smyth

I have some control spots that are not differentially expressed on the array.
I tried to use the modifyWeights as on page 21 of the Users Guide (6 Oct 2005).
However I realized that it doesn't work, cause I don't use weights in my
read.maimages. I used a image analysis program, not recognized by limma, so used
the following to read the data.

RG <- read.maimages(targets$FileName, columns=list(R="Rf",G="Gf",Rb="Rb",Gb="Gb"))

I tried to add "phony" column for weights (all values=1) in the intensity files,
so I can modify it, but have not been able to get it to work (yet).

Cecilia 

----- Original Message -----
From: "Gordon K Smyth" <smyth at wehi.EDU.AU>
To: "Cecilia McGregor" <cmcgre1 at lsu.edu>
Subject: [BioC] limma Normalization question
Date: Tue, 22 Nov 2005 07:54:59 +1100 (EST)

> 
> 
> > Date: Sun, 20 Nov 2005 17:09:09 -0600
> > From: "Cecilia McGregor" <cmcgre1 at lsu.edu>
> > Subject: [BioC] limma Normalization question
> > To: bioconductor at stat.math.ethz.ch
> >
> > Hi Everyone
> >
> > I've described my experiment in an earlier message, on Nov 16. 
> > After running the commands
> > (see end of message for commands) I looked at the 
> > plotMA(fit2)plot. I attached the plot.
> > It seems to me that at high A-values, the plot is going more in 
> > the direction of positive
> > M-values. I tested a few genes in the circled area (for high 
> > A-values)with Q-RT-PCR and it
> > confirms that these genes should have M-values around 0. The fact 
> > that more genes
> > are down-regulated, than up-regulated is expected from our 
> > knowledge of the experiment.
> > In the limma Users Guide (6 Oct 2005) p21, it says that loess 
> > does not assume equal number
> > of genes up- and down-regulated, so this should not be a problem. 
> > I have about 2500 unique
> > genes on the array.However it does say that most genes need to be 
> > not differentially
> > expressed. I estimate about 60% of genes are not differentially expressed.
> 
> This is really pushing the envelope.  I guess that loess can 
> probably handle up to around 20%
> differentially expressed, while robustspline can handle around 30%. 
>   With 40% you probably need
> some symmetry as well.
> 
> > In a message
> > on this board (Sept 3, 2003), Gordon Smyth suggests the use of 
> > the following command if
> > a large number of genes are differentially expressed.
> >
> > normalizeWithinArrays(RG, method="robustspline", robust="MM")
> >
> > However I get 18 of the following warning messages:
> > Warning messages:
> > 1: rlm failed to converge in 20 steps in: rlm.default(x, y, 
> > weights = w, method = method)
> 
> This is probably not a problem, so doesn't need a fix.
> 
> Another alternative would be to increase the robustifying 
> iterations used with "loess", say
> 
> > MA <- normalizeWithinArrays(RG, method="loess", iterations=10)
> 
> That is more robust than the default "loess" method.
> 
> Do you have any useful control probes on your arrays?
> 
> Best wishes
> Gordon
> 
> > I believe Paul Boutros suggested a fix for this in a message on 
> > this board(Aug 11 2004),
> > but since I don't know much about the subject, I'm not very eager 
> > to change the limma
> > code (I think that is what he suggest).
> >
> > How can I fix this problem I'm having at high A-Values? I get the 
> > correct results at low
> > A-values. I have analised this data also with limmaGUI (I just 
> > treated all replicates
> > as biological)and I get the similar results. I used several 
> > different within and between
> > normalization methods, but cannot get better results.
> >
> > I even tried normalizing with Genespring, altough I know it is 
> > not apropriate for loop
> > designs. With Genesspring I can get an acceptable MA plot if I 
> > use "Per array normalizations
> > to the 20th percentile". However I'm not very comfortable with 
> > this, since inspecting the
> > MA plots before normalization, I can see that I have an intensity 
> > dependent bias on all my
> > arrays.
> >
> > Any suggestions would be much appreciated.
> >
> > Code I Used: (I get the same results with the code suggested by 
> > Steen Krogsgaard on Nov 17)
> >
> > library(limma)
> > setwd("C:\\Program Files\\R_JSM\\rw2011\\RFlimma")
> > targets <- readTargets()
> > show(targets)
> > RG <- read.maimages(targets$FileName,
> > columns=list(R="Rf",G="Gf",Rb="Rb",Gb="Gb"))
> > RG$genes <- readGAL()
> > spottypes <-readSpotTypes()
> > RG$genes$Status <-controlStatus(spottypes, RG)
> > MA <- normalizeWithinArrays(RG, method="loess")
> > MA <- normalizeBetweenArrays(MA, method="Aquantile")
> > design <-modelMatrix(targets,ref="S1")
> > cor <- duplicateCorrelation(MA,design,ndups=3,spacing=3120)
> > cor$consensus.correlation
> > fit <-lmFit(MA,design,ndups=3,spacing=3120,correlation=0.3128828)
> > cont.matrix <- makeContrasts(fvss=(F1+F2+F3-S2-S3)/3, levels=design)
> > fit2 <- contrasts.fit(fit, cont.matrix)
> > fit2 <- eBayes(fit2)
> > plotMA(fit2)
> >
> > Thanks
> > Cecilia McGregor
> >
> > PhD Student
> > Sweetpotato Breeding and Genetics Lab
> > JC Miller Hall room 236
> > Louisiana State University
> > Baton Rouge
> > LA, 70803
> > USA
> >
> > Phone: (225) 578 2173



Cecilia McGregor

PhD Student
Sweetpotato Breeding and Genetics Lab
JC Miller Hall room 236
Louisiana State University
Baton Rouge
LA, 70803
USA

Phone: (225) 578 2173