[BioC] Limma, decideTests

[Ingunn Berget]

Hello mail-list

Experiment
This microarray experiment was conducted to study gene expression in 
bacteria at different growth conditions. There are 12 different conditions:

> unique(conditions)
 [1] "gr46"     "NaCl"     "Glycerol" "HCl"      "NaOH"     "CH3COOH"
 [7] "EtOH"     "gr15"     "PK"       "BC7ppm"   "BC9ppm"   "EtBr"

The one called PK is the positive control (normal conditions), whereas the 
other conditions are "stress conditions". A common reference, used for all 
arrays, is labelled with Cy3. On each array, one of the other conditions is 
labelled with Cy5. There are three biologial replicates of each 
array/condition, expect NaOH and HCl where something went wrong with one of 
replicates. Totally 34 arrays.

The goal is to find which genes are differentially expressed in PK and the 
other conditions.

ANALYSIS, copy of RCODE + some comments:
MAptl2 <- normalizeWithinArrays(modobj2)#modobj2: RGlist with raw data
D <- modelMatrix(modobj2$target,ref = "Ref")
GG <- MAptl2[keep,] #keep: logical index because want to use genes only, not 
controls
cor<- duplicateCorrelation(GG,design=D)
fit <- lmFit(GG,design=D,ndups=2,correlation = cor$consensus.correlation)
fit <- eBayes(fit)
topTable(fit,n=30,adjust="fdr")


The topTable command results in a list with genes having small p-values, 
have also tried with
topTable(fit,coef.=i,n=30,adjust="fdr") and different values of i 
(i=1,2,..,12)

> contrast.matrix <- 
> makeContrasts(gr15.PK=gr15-PK,gr46.PK=gr46-PK,gr46.15=gr46-gr15,NaCl.PK=NaCl-PK,EtBr.PK=EtBr-PK,
+ 
EtOH.PK=EtOH-PK,NaOH.PK=NaOH-PK,CH2COOH.PK=CH3COOH-PK,HCl.PK=HCl-PK,BC7.PK=BC7ppm-PK,
+ BC9.PK=BC9ppm-PK,BC9.BC7=BC9ppm-BC7ppm,Gly.PK=Glycerol-PK,levels=D)
>
> fit2 <- contrasts.fit(fit,contrast.matrix)
> fit2 <- eBayes(fit2)
>
> results <- decideTests(fit2,method="global")
> summary(results)

   gr15.PK gr46.PK gr46.15 NaCl.PK EtBr.PK EtOH.PK NaOH.PK CH2COOH.PK HCl.PK 
BC7.PK BC9.PK BC9.BC7 Gly.PK
-1
0
1
> table(results)
results
   -1     0     1
   13 78780    26
>


Questions:
Since the topTable command lead to a genelist with very low p-values, this 
means that there are differentially expressed genes in the data?
This is differentially expressed compared to the common reference?
Is this contrast matrix the correct for comparing gene expression in the 
stress conditions and PK
I don't understand why the summary(result) is "empty", does this means that 
a)  there's something wrong in the code, b)  there is not enough data for 
making this many contrasts (only 3 replicates?) c) no differential expressed 
genes??

If it is of help, here is a part of the target-slot in the original RG list 
(have not indluded all since this mail already is long)

   Slidenumber               FileNameCy3               FileNameCy5      Cy5 
Cy3         Name
1         1279     Cy3_1279_46gr_III.txt     Cy5_1279_46gr_III.txt     gr46 
Ref     46gr_III
2         1608      Cy3_1608_NaCl_II.txt      Cy5_1608_NaCl_II.txt     NaCl 
Ref      NaCl_II
3         1609  Cy3_1609_Glycerol_II.txt  Cy5_1609_Glycerol_II.txt Glycerol 
Ref  Glycerol_II
4         1610      Cy3_1610_HCl_III.txt      Cy5_1610_HCl_III.txt      HCl 
Ref      HCl_III
5         1612     Cy3_1612_NaOH_III.txt     Cy5_1612_NaOH_III.txt     NaOH 
Ref     NaOH_III
6         1613  Cy3_1613_CH3COOH_III.txt  Cy5_1613_CH3COOH_III.txt  CH3COOH 
Ref  CH3COOH_III
7         1625      Cy3_1625_46gr_II.txt      Cy5_1625_46gr_II.txt     gr46 
Ref      46gr_II

all 34 rows have Ref in the Cy3 column, in the Cy5 column according to the 
condition RNA is extracted from on each array

Thanks for any help in advance, and sorry for long mail

Best regards

Ingunn