[BioC] edgeR glm

[Gordon K Smyth]

Dear Shreyartha,

You code is actually testing for genes that are DE between f1 and p1.  By 
default, glmLRT tests for the last coefficient in the linear model and in 
your model this happens to represent f1-p1.  See help("glmLRT").

What you probably want is

   lrt <- glmLRT(d,fit,coef=c(2,3))

This will test whether f1=p1 and p2=p1 simultaneously, i.e., whether all 
three genotypes have the same expression.  This will find genes that are 
differentially expressed between any of the genotypes.  This is analogous 
to doing genewise oneway ANOVA tests for differences between the three 
genotypes.

Note that the above test does not guarantee that all three genotypes are 
different.  There is actually no statistical test that can do that. 
Rather it detects genes for which the genotypes are not all equal.

Unfortunately, the current release version of topTags() won't work 
properly when you test for several coefficients as above.  We'll fix 
this for the next Bioconductor release in a couple of weeks time.  In the 
meantime, you could use the following code:

   o <- order(lrt$table$p.value)
   ntags <- 10
   lrt$table[o[1:ntags],]

to see the most significant tags.

Best wishes
Gordon

---------------------------------------------
Professor Gordon K Smyth,
NHMRC Senior Research Fellow,
Bioinformatics Division,
Walter and Eliza Hall Institute of Medical Research,
1G Royal Parade, Parkville, Vic 3052, Australia.
Tel: (03) 9345 2326, Fax (03) 9347 0852,
smyth at wehi.edu.au
http://www.wehi.edu.au
http://www.statsci.org/smyth

On Tue, 5 Apr 2011, Shreyartha Mukherjee wrote:

> Hi Gondon,Bioconductor experts,
> 
> I have RNA-seq data for 3 genotypes(each genotype having 3 biological 
> replicates) and 30,000 genes. I was trying to find out differentially 
> expressed genes across all genotypes. Here is the code I am using.
>
> library(edgeR)
> y<-read.table('test_12339.txt'
> )
> d <- DGEList(y, group = rep(1:3, each = 3), lib.size = lib.sizes)
> d <- calcNormFactors(d)
> times <- rep(c("p1","p2","f1"),c(3,3,3))
> times <- factor(times,levels=c("p1","p2","f1"))
> design <- model.matrix(~factor(times))
> disp <- d2$tagwise.dispersion
> fit <- glmFit(d,design,dispersion=disp)
> lrt <- glmLRT(d,fit)
> topTags(lrt)
>
> Does this code provide me with tags differentially expressed across all
> conditions? (I am not looking for pairwise differential expression
> but across all conditions)
>
> Thanks,
> Shreyartha

______________________________________________________________________
The information in this email is confidential and intend...{{dropped:4}}