how to design a model matrix and test by edgeR dear all: my data is composed of 35 samples, have two factor time and treatment 1. The first question is i want to find DE genes (6h.treatment-0h.treatment)-(6h.control- 0h.control) how to design the model? please tell me if the design and test is right such is my coding. 2. the second question is should i use samples involved 6h and 0h to do glmFitting or use all of 35 samples to do? i think i should use all raw.data <- read.table("expression-table.txt",row.names=1) lib_size <- read.table("lib_size.txt"); lib_size <- unlist(lib_size) d <- DGEList(counts = raw.data, lib.size = lib_size) #normalization dge <- calcNormFactors(dge) treatment=factor(c(rep('control',6),rep('treated',24),rep('control',5) )) time=factor(c('0h','0h','0h','24h','24h','24h','0h','0h','0h','6h','6h ','6h','6h','12h','12h','12h','12h','18h','18h','18h','18h', '24h','24h','24h','36h','36h','36h','48h','48h','48h','6h ','12h','18h','36h','48h')) design <- model.matrix(~time*treatment) dge <- estimateGLMCommonDisp(dge, design) dge <- estimateGLMTagwiseDisp(dge, design) glmfit.dge <- glmFit(dge, design,dispersion=dge$common.dispersion) lrt.dge <- glmLRT(dge, glmfit.dge, coef=13) result <- topTags(lrt.dge, adjust.method="BH", sort.by="logFC") -- shan gao Room 231(Dr.Fei lab) Boyce Thompson Institute Cornell University Tower Road, Ithaca, NY 14853-1801 Office phone: 1-607-254-1267(day) Official email:sg839 at cornell.edu Facebook:http://www.facebook.com/profile.php?id=100001986532253

> Date: Wed, 25 Jul 2012 22:07:35 -0400 > From: wang peter <wng.peter at="" gmail.com=""> > To: bioconductor at r-project.org > Subject: [BioC] how to design a model matrix and test by edgeR > > how to design a model matrix and test by edgeR > dear all: > my data is composed of 35 samples, have two factor > time and treatment > > 1. The first question is > i want to find DE genes (6h.treatment-0h.treatment)-(6h.control- 0h.control) > how to design the model? please tell me if the design and test is right > such is my coding. As you can see for yourself, coef=13 does not relate to 0h: colnames(design)[13] [1] "time48h:treatmenttreated" On the other hand, coef=14 does: colnames(design)[14] [1] "time6h:treatmenttreated" So you could use either glmLRT(dge, glmfit.dge, coef=14) or glmLRT(dge, glmfit.dge, coef="time6h:treatmenttreated") to find the DE genes for the interaction you want. > 2. the second question is > should i use samples involved 6h and 0h to do glmFitting > or use all of 35 samples to do? > i think i should use all Yes, it is better to use all. There are some strange aspects to your code. You use common dispersion, although the edgeR User's Guide and our published papers make it clear that we recommend tagwise dispersions. You actually estimate tagwise dispersions, but then don't use them. Are you sure you need to sort genes by logFC? Best wishes Gordon > raw.data <- read.table("expression-table.txt",row.names=1) > lib_size <- read.table("lib_size.txt"); > lib_size <- unlist(lib_size) > d <- DGEList(counts = raw.data, lib.size = lib_size) > #normalization > dge <- calcNormFactors(dge) > > treatment=factor(c(rep('control',6),rep('treated',24),rep('control', 5))) > time=factor(c('0h','0h','0h','24h','24h','24h','0h','0h','0h','6h',' 6h','6h','6h','12h','12h','12h','12h','18h','18h','18h','18h', > '24h','24h','24h','36h','36h','36h','48h','48h','48h','6 h','12h','18h','36h','48h')) > design <- model.matrix(~time*treatment) > > dge <- estimateGLMCommonDisp(dge, design) > dge <- estimateGLMTagwiseDisp(dge, design) > glmfit.dge <- glmFit(dge, design,dispersion=dge$common.dispersion) > lrt.dge <- glmLRT(dge, glmfit.dge, coef=13) > result <- topTags(lrt.dge, adjust.method="BH", sort.by="logFC") > > -- > shan gao > Room 231(Dr.Fei lab) > Boyce Thompson Institute > Cornell University > Tower Road, Ithaca, NY 14853-1801 > Office phone: 1-607-254-1267(day) > Official email:sg839 at cornell.edu > Facebook:http://www.facebook.com/profile.php?id=100001986532253 ______________________________________________________________________ The information in this email is confidential and intend...{{dropped:4}}