[BioC] use of duplicateCorrelation in Limma with agilent one-color arrays

Tue Feb 21 04:19:39 CET 2012

Well, you could simplify to

Obj <- read.maimages(targets, source="agilent.median", green.only=TRUE)

and I recommend something more like offset=16 (any offset will give 
positive results).  Also coef=-1 must surely give an empty result.

Best wishes
Gordon

---------------------------------------------
Professor Gordon K Smyth,
Bioinformatics Division,
Walter and Eliza Hall Institute of Medical Research,
1G Royal Parade, Parkville, Vic 3052, Australia.
smyth at wehi.edu.au
http://www.wehi.edu.au
http://www.statsci.org/smyth

On Mon, 20 Feb 2012, Jabez Wilson wrote:

> Dear Gordon, thank you for your comments. I am new to single color arrays, having used limma extensively with two-color in the past. I admit I have borrowed heavily with my coding from the mattick lab wiki(http://matticklab.com/index.php?title=Single_channel_analysis_of_Agilent_microarray_data_with_Limma), which you may have a view on.
> I would be grateful if you could comment therefore on my amended script? (which with my data gives a corfit$consensus of ~0.4)
> (also assuming the design matrix is correctly constructed)
> Obj <- read.maimages(targets, source="agilent", green.only=TRUE, columns=list(G="gMedianSignal", Gb="gBGMedianSignal"))
> Obj.corr <- backgroundCorrect (Obj, method="normexp", offset=1)
> E <- normalizeBetweenArrays(Obj.corr, method="quantile")
> #E.avg <= avereps(E, ID=E$genes$ProbeName)
> then the analysis similar to before using E instead of E.avg
> biolrep <- c(1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8,9.9,10,10,11,11,12,12,13,13,14,14,15,15,16,16)
> corfit <- duplicateCorrelation(E, design, ndups=1, block=biolrep)
> fit <- lmFit(E, design, block=biolrep, cor=corfit$consensus)
> cont.matrix <- makeContrasts("Treat1-WT", levels=design)
> fit2 <- contrasts.fit(fit, cont.matrix)
> fit2 <- eBayes(fit2)
> topTable(fit2, coef=-1)
>
> Jabez
>
>
>
>
> --- On Sun, 19/2/12, Gordon K Smyth <smyth at wehi.EDU.AU> wrote:
>
>> From: Gordon K Smyth <smyth at wehi.EDU.AU>
>> Subject: use of duplicateCorrelation in Limma with agilent one-color arrays
>> To: "Jabez Wilson" <jabezwuk at yahoo.co.uk>
>> Cc: "Bioconductor mailing list" <bioconductor at r-project.org>
>> Date: Sunday, 19 February, 2012, 10:20
>> Dear Jabez,
>>
>> There a number of strange things about your code, which seem
>> to be to do with trying to work around storing single color
>> data in a 2-color data object.  Could you please use
>> read.maimages() with green.only=TRUE, so that the data is
>> read into a single color data object.  None of the
>> work-around will then be necessary.
>>
>> You might also try computing the duplicate correlation
>> without averaging replicate probes.
>>
>> Best wishes
>> Gordon
>>
>>> Date: Fri, 17 Feb 2012 12:14:23 +0000 (GMT)
>>> From: Jabez Wilson <jabezwuk at yahoo.co.uk>
>>> To: bioconductor at r-project.org
>>> Subject: [BioC] use of duplicateCorrelation in Limma
>> with agilent
>>>     one-color    arrays
>>> Message-ID:
>>>     <1329480863.2527.YahooMailClassic at web132503.mail.ird.yahoo.com>
>>> Content-Type: text/plain; charset=iso-8859-1
>>>
>>> Hi, everyone, I am using limma to analyse an agilent
>> one-color array experiment, and have run into difficulties
>> with duplicateCorrelation.
>>> My experiment is as follows: single color agilent
>> arrays, 4 WT samples, and 3 samples of each of 4 treatment
>> (treatments 1-4). I also have technical replicates
>> (replicated once) for each sample. There are therefore 32
>> files. The targets file looks like this:
>>> ?
>>> ?? SampleNumber????????????????? FileName Condition
>> Notes
>>> 1???????????? 1? RH_02_1_77_Oct11_1_1.txt??? Treat1??
>> 1.1
>>> 2???????????? 2? RH_02_1_77_Oct11_2_1.txt??? Treat1??
>> 1.2
>>> 3???????????? 3? RH_04_1_77_Oct11_1_1.txt??? Treat1??
>> 2.1
>>> 4???????????? 4? RH_07_1_77_Oct11_1_1.txt??? Treat1??
>> 2.2
>>> 5???????????? 5? RH_04_1_77_Oct11_1_2.txt??? Treat1??
>> 3.1
>>> 6???????????? 6? RH_07_1_77_Oct11_1_2.txt??? Treat1??
>> 3.2
>>> 7???????????? 7? RH_04_1_77_Oct11_1_3.txt??? Treat2??
>> 4.1
>>> 8???????????? 8? RH_07_1_77_Oct11_1_3.txt??? Treat2??
>> 4.2
>>> 9???????????? 9? RH_04_1_77_Oct11_1_4.txt??? Treat2??
>> 5.1
>>> 10?????????? 10? RH_07_1_77_Oct11_1_4.txt??? Treat2??
>> 5.2
>>> 11?????????? 11? RH_04_1_77_Oct11_2_1.txt??? Treat2??
>> 6.1
>>> 12?????????? 12? RH_07_1_77_Oct11_2_1.txt??? Treat2??
>> 6.2
>>> 13?????????? 13 US0_05_1_77_Oct11_1_1.txt??? Treat3??
>> 7.1
>>> 14?????????? 14? RH_01_1_77_Oct11_1_1.txt??? Treat3??
>> 7.2
>>> 15?????????? 15 US0_05_1_77_Oct11_1_2.txt??? Treat3??
>> 8.1
>>> 16?????????? 16? RH_01_1_77_Oct11_1_4.txt??? Treat3??
>> 8.2
>>> 17?????????? 17? RH_02_1_77_Oct11_1_2.txt??? Treat3??
>> 9.1
>>> 18?????????? 18? RH_02_1_77_Oct11_2_2.txt??? Treat3??
>> 9.2
>>> 19?????????? 19 US0_05_1_77_Oct11_1_3.txt??? Treat4?
>> 10.1
>>> 20?????????? 20? RH_01_1_77_Oct11_1_2.txt??? Treat4?
>> 10.2
>>> 21?????????? 21 US0_05_1_77_Oct11_1_4.txt??? Treat4?
>> 11.1
>>> 22?????????? 22? RH_01_1_77_Oct11_1_3.txt??? Treat4?
>> 11.2
>>> 23?????????? 23? RH_04_1_77_Oct11_2_2.txt??? Treat4?
>> 12.1
>>> 24?????????? 24? RH_07_1_77_Oct11_2_2.txt??? Treat4?
>> 12.2
>>> 25?????????? 25 US0_05_1_77_Oct11_2_1.txt??????? WT?
>> 13.1
>>> 26?????????? 26? RH_01_1_77_Oct11_2_1.txt??????? WT?
>> 13.2
>>> 27?????????? 27 US0_05_1_77_Oct11_2_2.txt??????? WT?
>> 14.1
>>> 28?????????? 28? RH_01_1_77_Oct11_2_2.txt??????? WT?
>> 14.2
>>> 29?????????? 29 US0_05_1_77_Oct11_2_3.txt??????? WT?
>> 15.1
>>> 30?????????? 30? RH_01_1_77_Oct11_2_3.txt??????? WT?
>> 15.2
>>> 31?????????? 31 US0_05_1_77_Oct11_2_4.txt??????? WT?
>> 16.1
>>> 32?????????? 32? RH_01_1_77_Oct11_2_4.txt??????? WT?
>> 16.2
>>>
>>> I run the following commands to process the data and
>> create the design:
>>> RG <- read.maimages(targets, columns = list(G =
>> "gMedianSignal", Gb = "gBGMedianSignal", R =
>> "gProcessedSignal",Rb = "gIsPosAndSignif"), annotation =
>> c("Row", "Col","FeatureNum","ControlType","ProbeName"))
>>> RG <- backgroundCorrect(RG, method="normexp",
>> offset=1)
>>> E <- normalizeBetweenArrays(RG, method="Aquantile")
>>> E.avg <- avereps(E, ID=E$genes$ProbeName)
>>> f <- factor(targets$Condition, levels =
>> unique(targets$Condition))
>>> design <- model.matrix(~0 + f)
>>> colnames(design) <- levels(f)
>>> ?
>>> The problem arises when I do the duplicateCorrelation.
>>> biolrep <-
>> c(1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8,9,9,10,10,11,11,12,12,13,13,14,14,15,15,16,16)
>>> corfit <- duplicateCorrelation(E.avg, design,
>> ndups=1,block=biolrep)
>>> fit <- lmFit(E.avg$A, design, block=biolrep,
>> cor=corfit$consensus.correlation)
>>> contrast.matrix <-
>> makeContrasts("Treat1-WT",levels=design)????????????????????????????????????????????????????????????fit2
>> <- contrasts.fit(fit, contrast.matrix)
>>> fit2 <- eBayes(fit2)
>>> topTable(fit2, adjust="BH", coef="Treat1-WT",
>> genelist=E.avg$genes, number=10)
>>>
>>> Whereas I would expect the corfit$consensus.correlation
>> to be generally very positive, I get the value
>> 0.01385223.?Does anyone have any suggestions? Any help would
>> be appreciated
>>> ?
>>> Jabez
>>>
>>
>> ______________________________________________________________________
>> The information in this email is confidential and intended
>> solely for the addressee.
>> You must not disclose, forward, print or use it without the
>> permission of the sender.
>> ______________________________________________________________________
>>
>

______________________________________________________________________
The information in this email is confidential and intend...{{dropped:5}}