dear all: my data is composed of 35 samples, have two factor time and treatment i want to find DE genes between 6h and 0h due to the treatment effect,removing the time factor effect. how can i use glamLRT to do test 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) dge <- d[rowSums(d$counts) >= length(lib_size)/2,] #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=2)?????????????????????????????? how can i use glamLRT to do test -- 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 -- 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

Hi Shan Gao, On 7/18/2012 10:34 AM, wang peter wrote: > dear all: > > my data is composed of 35 samples, have two factor > time and treatment > > i want to find DE genes between 6h and 0h due to the treatment > effect,removing the time factor effect. I am not sure what you mean. How can you find DE genes between two times after removing the time factor effect? The default of glmLRT(dge, glmfit.dge) will give you the treatment vs control comparison, adjusted for time. > > how can i use glamLRT to do test > > 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) > dge<- d[rowSums(d$counts)>= length(lib_size)/2,] > > #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=2)?????????????????????????????? This will give you the 12hr - 0 hr comparison, adjusted for treatment status. Best, Jim > > > > how can i use glamLRT to do test > > -- > 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 > > -- James W. MacDonald, M.S. Biostatistician University of Washington Environmental and Occupational Health Sciences 4225 Roosevelt Way NE, # 100 Seattle WA 98105-6099