Hello everyone, I'm trying to compare limma vs edgeR in finding tissue specific genes from my RNAseq samples. For edgeR I've followed directions from a previous post and have managed to get good results. The thread is named "edgeR: finding tissue specific genes" and here is the link: http://permalink.gmane.org/gmane.science.biology.informatics.conductor /47950 However, I'm having a difficult time creating the contrast for the last tissue in limma. For the last tissue, the contrast is a vector with 18 elements which lmFit does not like. First of all, here is a snippet of code of what I did in edgeR contrasts(tiss_groups) <- contr.sum(tiss_groups) design <- model.matrix(~tiss_groups) y <- estimateGLMCommonDisp(y, design) y <- estimateGLMTrendedDisp(y, design) fit <- glmFit(y, design) ## Get enriched genes for 1st tissue: ql <- glmQLFTest(fit, coef=2) top1 <- topTags(ql) ## Getting tissue 2 - 17 is the same as above but with different coef values ### #Get enrichment for the last tissue cont <- rep(-1,18) cont[1] <- 0 ql <- glmQLFTest(fit, contrast=cont) top18 <- topTags(de) Here is what I'm doing in limma: contrasts(tiss_groups) <- contr.sum(tiss_groups) design <- model.matrix(~tiss_groups) v <- voom(y, design, plot=T) fit <- lmFit(v,design) fit <- eBayes(fit) ##Get top genes in tissue j = 1:17 top1 <- topTable(fit, coef=(j+1), adjust="BH", sort="p") ## How to get top genes the last tissue? Here is what I did: cont <- rep(-1, 18) cont[1] <- 0 fit <- lmFit(v,design=cont) ##<------- this fails with error "(subscript) logical subscript too long" fit <- eBayes(fit) top18 <- topTable(fit, n=nrow(enrich), sort="p") Thanks for your help, Jaaved Mohammed Graduate student of Computational Biology Cornell University

Dear Jaaved, limma and edgeR work very similarly. If you use coef=2 in edgeR, then that's just what you use in limma also. You only have 17 tissues, so obviously trying to refer to coefficient number 18 in limma will cause an error. There is no tissue 18. Best wishes Gordon > Date: Tue, 25 Jun 2013 13:18:31 -0400 > From: Jaaved Mohammed <jm889 at="" cornell.edu=""> > To: <bioconductor at="" stat.math.ethz.ch=""> > Subject: [BioC] limma: tissue specific genes using voom > > Hello everyone, > > I'm trying to compare limma vs edgeR in finding tissue specific genes > from my RNAseq samples. For edgeR I've followed directions from a > previous post and have managed to get good results. The thread is named > "edgeR: finding tissue specific genes" and here is the link: > http://permalink.gmane.org/gmane.science.biology.informatics.conduct or/47950 > > However, I'm having a difficult time creating the contrast for the last > tissue in limma. For the last tissue, the contrast is a vector with 18 > elements which lmFit does not like. > > First of all, here is a snippet of code of what I did in edgeR > > contrasts(tiss_groups) <- contr.sum(tiss_groups) > design <- model.matrix(~tiss_groups) > y <- estimateGLMCommonDisp(y, design) > y <- estimateGLMTrendedDisp(y, design) > fit <- glmFit(y, design) > > ## Get enriched genes for 1st tissue: > ql <- glmQLFTest(fit, coef=2) > top1 <- topTags(ql) > > ## Getting tissue 2 - 17 is the same as above but with different coef > values ### > > #Get enrichment for the last tissue > cont <- rep(-1,18) > cont[1] <- 0 > ql <- glmQLFTest(fit, contrast=cont) > top18 <- topTags(de) > > > > > Here is what I'm doing in limma: > > contrasts(tiss_groups) <- contr.sum(tiss_groups) > design <- model.matrix(~tiss_groups) > > v <- voom(y, design, plot=T) > fit <- lmFit(v,design) > fit <- eBayes(fit) > > ##Get top genes in tissue j = 1:17 > top1 <- topTable(fit, coef=(j+1), adjust="BH", sort="p") > > ## How to get top genes the last tissue? Here is what I did: > cont <- rep(-1, 18) > cont[1] <- 0 > fit <- lmFit(v,design=cont) ##<------- this fails with error > "(subscript) logical subscript too long" > fit <- eBayes(fit) > top18 <- topTable(fit, n=nrow(enrich), sort="p") > > > > Thanks for your help, > Jaaved Mohammed > Graduate student of Computational Biology > Cornell University ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}