Search
Question: enquiry about the limma script
0
gravatar for yong li
5.7 years ago by
yong li10
yong li10 wrote:
Dear Sir/Madam, I meet some problems when I use the limma package to analysis the Agilent single channel microarray. I got 8 samples (with 8*60K array) results as following: Samples: (T1, T2, T3) =Treat1 Samples: (T4, T5, T6)=Treat2 Samples: (C1, C2)=Control my targets.txt contain: file_name conditon x1.txt T1 x2.txt T2 x3.txt T3 x4.txt T4 x5.txt T5 x6.txt T6 x7.txt C1 x8.txt C8 I want to get the different expressed genes under 3 conditions: Treat1 vs Control, Treat2 vs Control, Treat2 vs Treat1 with the following format: *GROUP1* *GROUP2* *CONTROL* *GROUP1-vs-CONTROL* *GROUP2-vs-CONTROL* *GROUP1-vs-GROUP2* *Probeid* *Accession* *T1* *T2* *T3* *T4* *T5* *T6* *C1* *C2* *avg-group1* *avg-group2* *avg-control* *FC* *UP/DOWN* *P value* *FC* *UP/DOWN* *P value* *FC* *UP/DOWN* *P value* CUST_PI5 AB91 3 1 1 1 2 4 6 6 2 2.5 6 6 DOWN 0.01 2 DOWN 0.02 1 DOWN 0.2 My R script as follow: library(limma) targets <- readTargets("targets.txt") x <- read.maimages(targets[,"file_name"], source="agilent",green.only=TRUE) y <- backgroundCorrect(x, method="normexp") y <- normalizeBetweenArrays(y, method="quantile") y.ave <- avereps(y, ID=y$genes$ProbeName) f <- factor(targets$condition, levels = unique(targets$condition)) design <- model.matrix(~0 + f) colnames(design) <- levels(f) fit <- lmFit(y.ave, design) contrast.matrix <- makeContrasts("Treat1-C", "Treat2-C","Treat2-Treat1",levels=design) fit2 <- contrasts.fit(fit, contrast.matrix) fit2 <- eBayes(fit2) write.table(topTable(fit2,adjust="BH",number="44000"), file="out_file.txt",sep="\t") I got the results as following: *ProbeName* *SystematicName* *Treat1.C* *Treat2.C* *Treat1.Treat2* *AveExpr* *F* *P.Value* *adj.P.Val* CUST_248 EV478 -0.54306 -8.90559 8.362525 14.51936 2562.081 1.97E-10 6.12E-06 The problem is that I can’t get the normalized data, every group average data, FC with up/down label in a file just like the previous format. I write you for help and look forward to receiving your reply. Thanks a lot. Yong Li Institute of Plant Physiology and Ecology, SIBS, CAS [[alternative HTML version deleted]]
ADD COMMENTlink modified 5.7 years ago by Gordon Smyth35k • written 5.7 years ago by yong li10
0
gravatar for Gordon Smyth
5.7 years ago by
Gordon Smyth35k
Walter and Eliza Hall Institute of Medical Research, Melbourne, Australia
Gordon Smyth35k wrote:
Dear Yong Li, The formatted table you have asked for is not one produced automatically by limma. The closest is that produced by write.fit(). If you are trying to reproduce a previous analysis, you should contact the person who produced that the previous analysis. Best wishes Gordon > Date: Sun, 24 Mar 2013 20:25:20 +0800 > From: yong li <yonglinet at="" gmail.com=""> > To: <bioconductor at="" stat.math.ethz.ch=""> > Subject: [BioC] enquiry about the limma script > > Dear Sir/Madam, > > > > I meet some problems when I use the limma package to analysis the Agilent > single channel microarray. I got 8 samples (with 8*60K array) results as > following: > > Samples: (T1, T2, T3) =Treat1 > > Samples: (T4, T5, T6)=Treat2 > Samples: (C1, C2)=Control > my targets.txt contain: > file_name conditon > x1.txt T1 > x2.txt T2 > x3.txt T3 > x4.txt T4 > x5.txt T5 > x6.txt T6 > x7.txt C1 > x8.txt C8 > > > I want to get the different expressed genes under 3 conditions: Treat1 vs > Control, Treat2 vs Control, Treat2 vs Treat1 with the following format: > > *GROUP1* > > *GROUP2* > > *CONTROL* > > *GROUP1-vs-CONTROL* > > *GROUP2-vs-CONTROL* > > *GROUP1-vs-GROUP2* > > *Probeid* > > *Accession* > > *T1* > > *T2* > > *T3* > > *T4* > > *T5* > > *T6* > > *C1* > > *C2* > > *avg-group1* > > *avg-group2* > > *avg-control* > > *FC* > > *UP/DOWN* > > *P value* > > *FC* > > *UP/DOWN* > > *P value* > > *FC* > > *UP/DOWN* > > *P value* > > CUST_PI5 > > AB91 > > 3 > > 1 > > 1 > > 1 > > 2 > > 4 > > 6 > > 6 > > 2 > > 2.5 > > 6 > > 6 > > DOWN > > 0.01 > > 2 > > DOWN > > 0.02 > > 1 > > DOWN > > 0.2 > > > > My R script as follow: > > library(limma) > > targets <- readTargets("targets.txt") > > x <- read.maimages(targets[,"file_name"], source="agilent",green.only=TRUE) > > y <- backgroundCorrect(x, method="normexp") > > y <- normalizeBetweenArrays(y, method="quantile") > > y.ave <- avereps(y, ID=y$genes$ProbeName) > > f <- factor(targets$condition, levels = unique(targets$condition)) > > design <- model.matrix(~0 + f) > > colnames(design) <- levels(f) > > > > fit <- lmFit(y.ave, design) > > contrast.matrix <- makeContrasts("Treat1-C", > "Treat2-C","Treat2-Treat1",levels=design) > > fit2 <- contrasts.fit(fit, contrast.matrix) > > fit2 <- eBayes(fit2) > > write.table(topTable(fit2,adjust="BH",number="44000"), > file="out_file.txt",sep="\t") > > > > I got the results as following: > > *ProbeName* > > *SystematicName* > > *Treat1.C* > > *Treat2.C* > > *Treat1.Treat2* > > *AveExpr* > > *F* > > *P.Value* > > *adj.P.Val* > > CUST_248 > > EV478 > > -0.54306 > > -8.90559 > > 8.362525 > > 14.51936 > > 2562.081 > > 1.97E-10 > > 6.12E-06 > > > The problem is that I can?t get the normalized data, every group average > data, FC with up/down label in a file just like the previous format. I > write you for help and look forward to receiving your reply. > > Thanks a lot. > > > > Yong Li > > Institute of Plant Physiology and Ecology, SIBS, CAS > ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:5}}
ADD COMMENTlink written 5.7 years ago by Gordon Smyth35k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 365 users visited in the last hour