Dear all,

I'm analysing data single color data (only intensities at F532) using limma.

I have 14 arrays in one group and 5 in the other group. I want to compare them 
so my targets file is as follows
FileName Mal Sain
M2.gpr 1 0
M3.gpr 1 0
M4.gpr 1 0
M5.gpr 1 0
M6.gpr 1 0
M7.gpr 1 0
M8.gpr 1 0
M9.gpr 1 0
M10.gpr 1 0
M11.gpr 1 0
M12.gpr 1 0
M13.gpr 1 0
M14.gpr 1 0
M15.gpr 1 0
S1.gpr 0 1
S2.gpr 0 1
S3.gpr 0 1
S4.gpr 0 1
S5.gpr 0 1

I use some correction and normalisation methods as recommanded in limma user 
guide. Afterthat, I make my design and as I have 2 technical replicate on each 
array, I wanted to compute duplicatecorrelation coefficient with spacing=500 ( 
ref. .GAL file as follows : 

ATF 1
52 5
Type=GenePix ArrayList V1.0 
BlockCount=48 
BlockType=0 
URL=http://genome-www4.stanford.edu/cgi-bin/SGD/locus.pl?locus=[ID] 
"Block1=2645, 2295, 100, 12, 300, 8, 300" 
"Block2=7145, 2295, 100, 12, 300, 8, 300" 
"Block3=11555, 2235, 100, 12, 300, 8, 300" 
"Block4=16055, 2235, 100, 12, 300, 8, 300" 

The  corfit$consensus is not good (0,36). Do I still use this parameters for 
lmfit???

Anyway, the more strange is with or without this parameters I have p-adjust 
value very very low. 

first line of my result :
Block Row Column ID Name X.Intercept. factor.pData.population.sain AveExpr F P.Value adj.P.Val

39 2 10 42554 hsa-miR-923 13.4125293049973 -0.0517387609947947 13.3989138415776 25264.4487822867 1.20563134868457e-63 1.40456052121752e-60


Is anybody help me and what is wrong.

As you can see in UNFILTERED section, I have less than 50% genes for the 
analysis. Could it be the explanation of these strange result??

Thanks in advance

My srcipt:
library(limma)
myFilter = function(X) {
H_threshold=2
okFLAG = X$Flags > -49;
okm1 =abs(X[,"F532 Median"]-X[,"F532 Mean"])  
okm2 = 0.5*(X[,"F532 Median"]+X[,"F532 Mean"])
okH = ((okm1)/okm2) <H_threshold
as.numeric(okFLAG & okH)}

targets=readTargets(file="targets.txt", path=NULL, sep="\t")

E =read.maimages(targets,source="genepix",wt.fun=myFilter,columns=list(E="F532 
Mean",Eb="B532 Mean",Names="Name"))

#Taux de filtrage
unFiltered = apply(E$weights,MARGIN=2,FUN=mean)
round(unFiltered,2)  

#  M2   M3   M4   M5   M6   M7   M8   M9  M10  M11  M12  M13  M14  M15   S1   S2 

#0.41 0.43 0.34 0.29 0.51 0.35 0.48 0.35 0.29 0.56 0.38 0.56 0.61 0.52 0.48 0.34 

#  S3   S4   S5 
#0.36 0.46 0.47

#Boite a moustache avant correction
boxplot(as.data.frame(E$E),main="Mean intensities") 

#Correction background - subtract
Esub <- backgroundCorrect(E, method="subtract",offset=0)
boxplot(as.data.frame(Esub$E), main="Mean intensities - subtract")
#Correction background - minimum
Emin <- backgroundCorrect(E, method="minimum",offset=1)
boxplot(as.data.frame(Emin$E), main="Mean intensities - minimum")
#Correction background - normexp
Enorm <- backgroundCorrect(E, method="normexp",offset=1)
boxplot(as.data.frame(Enorm$E), main="Mean intensities - normexp")
#Normalisation quantile
NormE <- normalizeBetweenArrays(Enorm,method="quantile")
boxplot(as.data.frame(NormE$E), main="Normalized intensities")
pData <- data.frame(population = c('mal','mal', 'mal', 
'mal','mal','mal','mal','mal','mal','mal','mal','mal','mal','mal','sain','sain','sain','sain','sain'))


# create design matrix
design <- model.matrix(~factor(pData$population))

# In my .gpr files all miRNAs contain the string "miR" in their name
# so the grep function can be used to extract all of these, removing
# all control signals and printing buffers etc.
# You need to check your .gpr files to find which signals you should extract.
miRs <- c(grep("-miR-", NormE$genes$Name), grep("-let-", NormE$genes$Name))
NormE.final <- NormE[miRs, ]
NormE.final <-NormE.final[order(NormE.final$genes$Name), ]
NormE.final$genes$Name
# calculate duplicate correlation between the 2 replicates on each array
corfit <- duplicateCorrelation(NormE.final, design, ndups=2,spacing=500,)
corfit$consensus
# show a boxplot of the correlations
boxplot(tanh(corfit$atanh.correlations))

# fit the linear model, including info on duplicates and correlation
fit <- lmFit(NormE.final, design, ndups=2,  correlation=corfit$consensus)
#fit <- lmFit(NormE.final, design)
# calculate values using ebayes
ebayes <- eBayes(fit)

# output a list of top differnetially expressed genes...
result=topTable(ebayes,number=4608, p=0.05,adjust = "BH")
genes=result$Status=="genes"
write.table(result,"gene.txt",quote=FALSE,sep="\t",row.names=FALSE, 
col.names=TRUE)


mikel
	[[alternative HTML version deleted]]

