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/