[BioC] Design in factDesign

[Jordi Altirriba Gutiérrez]

Hi all!
I’ve been using RMA and LIMMA to analyse my data and I am currently trying 
to analyse it with the package factDesign. My design is a 2x2 factorial 
design with 4 groups: diabetic treated, diabetic untreated, health treated 
and health untreated with 3 biological replicates in each group. I want to 
know what genes are differentially expressed due to diabetes,  to the 
treatment and to the combination of both (diabetes + treatment).
My phenoData is:
>pData(eset)
         DIABETES    TREATMENT
DNT1     TRUE       FALSE
DNT2     TRUE       FALSE
DNT3     TRUE       FALSE
DT1        TRUE       TRUE
DT2        TRUE       TRUE
DT3        TRUE       TRUE
SNT1      FALSE     FALSE
SNT2      FALSE     FALSE
SNT3      FALSE     FALSE
ST1         FALSE     TRUE
ST2         FALSE     TRUE
ST3         FALSE     TRUE

Are these commands correct to get the results listed below? Where are the 
errors?
>lm.full<-function(y) lm(y ~ DIABETES + TREATMENT + DIABETES * TREATMENT)
>lm.diabetes<-function(y) lm(y ~ DIABETES)
>lm.treatment<-function(y) lm(y ~ TREATMENT)
>lm.diabetestreatment<-function(y) lm(y ~ DIABETES + TREATMENT)
>lm.f<-esApply(eset, 1, lm.full)
>lm.d<-esApply(eset, 1, lm.diabetes)
>lm.t<-esApply(eset, 1, lm.treatment)
>lm.dt<-esApply(eset, 1, lm.diabetestreatment)

## To get the genes characteristics of the treatment:
>Fpvals<-rep(0, length(lm.f))
>for (i in 1:length(lm.f)) {Fpvals[i]<-anova(lm.d[[i]], lm.f[[i]])$P[2]}
>Fsub<-which(Fpvals<0.01)
>eset.Fsub<-eset[Fsub]
>lm.f.Fsub<-lm.f[Fsub]
>betaNames<-names(lm.f[[1]] [["coefficients"]])
>lambda<-par2lambda(betaNames, c("TREATMENTTRUE"), c(1)) ## I get the same 
>genes if I write : > lambda2 <- par2lambda (betaNames, 
>list(c("TREATMENTTRUE" , "DIABETESTRUE:TREATMENTTRUE")),list( c(1,1)))
>mainTR<-function(x) contrastTest(x,lambda,p=0.1) [[1]]
>mainTRgenes<-sapply(lm.f.Fsub, FUN=mainES)

## To get the genes characteristics of the diabetes:
>for (i in 1:length(lm.f)) {Fpvals[i]<-anova(lm.t[[i]], lm.f[[i]])$P[2]}
>Fsub<-which(Fpvals<0.01)
>eset.Fsub<-eset[Fsub]
>lm.f.Fsub<-lm.f[Fsub]
>betaNames<-names(lm.f[[1]] [["coefficients"]])
>lambda<-par2lambda(betaNames, c("DIABETESTRUE"), c(1)) ## I get also the 
>same genes if I consider the intersection DIABETESTRUE:TREATMENTTRUE.
>mainDI<-function(x) contrastTest(x,lambda,p=0.1) [[1]]
>mainDIgenes<-sapply(lm.f.Fsub, FUN=mainES)

## To get the genes characteristics of the diabetes+treatment:
>for (i in 1:length(lm.f)) {Fpvals[i]<-anova(lm.dt[[i]], lm.f[[i]])$P[2]}
>Fsub<-which(Fpvals<0.01)
>eset.Fsub<-eset[Fsub]
>lm.f.Fsub<-lm.f[Fsub]
>  betaNames<-names(lm.f[[1]] [["coefficients"]])
>lambda<-par2lambda(betaNames, c("DIABETESTRUE:TREATMENTTRUE"), c(1))
>mainDT<-function(x) contrastTest(x,lambda,p=0.1) [[1]]
>mainDTgenes<-sapply(lm.f.Fsub, FUN=mainES) ## I don’t get any “fail to 
>reject” gene.

When I get the “rejected” and the “failed to reject” genes, can I classify 
them by their Fvalues? How?

Thank you very much for your comments and help.
Yours sincerely,

Jordi Altirriba
IDIBAPS-Hospital Clinic (Barcelona, Spain)

_________________________________________________________________
Déjanos tu CV y recibe ofertas de trabajo en tu buzón. Multiplica tus 
oportunidades con MSN Empleo. http://www.msn.es/Empleo/