[BioC] use of duplicateCorrelation in Limma with agilent one-color arrays
Jabez Wilson
jabezwuk at yahoo.co.uk
Fri Feb 17 13:14:23 CET 2012
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
More information about the Bioconductor
mailing list