Question: RNA-seq Anova on three groups
gravatar for BM
5.3 years ago by
United Kingdom
BM60 wrote:


I would like some advice on RNA-seq analysis. I have three groups (control, mutant 1 and mutant 2). we would like perform a one-way anova/GLM and to analyse all three groups pair wise. So analyse all three groups to test for significance, then esablish if any gene is changed bewteenany of the groups: control vs mutant 1; control vs mutant 2; mutant 1 vs mutant 2. However, we are running into problems.

Is it possible to do this in DESeq2 or EdgeR, or should iuse Genespring or Partek?


ancodev edger rna-seq deseq2 anova • 6.2k views
ADD COMMENTlink modified 5.3 years ago by Michael Love2.1k • written 5.3 years ago by BM60

Thanks. But what is the formula for doing this in DESeq2. So will Deseq2 analyse all three groups at once or only two (which I thought was the case). And also if all three are analysed in ANOVA/ANODEV, are pairwise compaisions performed an p-values listed for all pairwise comparisons or an overall comparision?

ADD REPLYlink written 5.3 years ago by BM60

Just read through the DESeq2 vignette. You would normally use a design like ~group and just specify each pairwise comparison in a different contrast. You'll then get 3 separate sets of adjusted p-values, just as you would with an ANOVA. You could also do an overall comparison if you wanted to, for which you would use a likelihood ratio test with a reduced model of ~1 (again, just read through the vignette for the exact commands).


ADD REPLYlink written 5.3 years ago by Devon Ryan97k

The data has been aligned using TopHat and counts performed using Dexseq.
I am trying to use Deseq2 for differential expression, but my experience of R  is limited.

1. I am having difficulty constructing the metadata on each of the samples (columns of the count table). I can load the sample info, but I can’t assign it to coldata.
So I have used this instead:

> counts = read.delim("dexseqcount.txt", header=T, row.names=1)
> sample <- read.delim("~/sample.txt")
> <- data.frame(genotype=as.factor(c("WT", "WT", "WT", "WT", "Mut1", "Mut1", "Mut1", "Mut2", "Mut2", "Mut2", "Mut2")))
> <- DESeqDataSetFromMatrix(countData=counts,,design= ~ genotype)
> dds<-DESeq(
> res <- results(dds)

2. Also if I want to a pairwise comparison between all three groups (WT, Mut1, Mut2), as the above is a comparison  between WT and Mut2 would this be ok:

> res2 <-results(dds, contrast=c("genotype","WT","Mut1"))
> res3 <-results(dds, contrast=c("genotype","Mut1","Mut2"))

3. I presume it is correct doing a pairwise comparison between all the three groups and looking to see if a gene is differentially expressed between WT v Mut2 and also in Mut1 vs Mut2 but not WT and Mut1. This will be done “by eye” by suing the three lists generated and seeing if there is any overlap

4. Final question, is it correct to do the above by doing the default test in Deseq2 i.e. Wald or should I use the alternative likelihood ratio test ?


ADD REPLYlink written 5.3 years ago by BM60

please do not cross post

ADD REPLYlink written 5.3 years ago by Michael Love2.1k
gravatar for Devon Ryan
5.3 years ago by
Devon Ryan97k
Freiburg, Germany
Devon Ryan97k wrote:

You can certainly do this in either DESeq2 or edgeR. You'll just specify the 3 contrasts, one for each comparison.

ADD COMMENTlink written 5.3 years ago by Devon Ryan97k
gravatar for Steven Lakin
5.3 years ago by
Steven Lakin1.5k
Fort Collins, CO, USA
Steven Lakin1.5k wrote:

Post-hoc testing is required once you have constructed your GLM.  The GLM will only estimate the dispersions amongst your variables; in order to get p-values, you have to perform testing based on the model.  Probably the easiest thing to do is to take a look at the edgeR vignette under the section 2.9 More complex experiments (glm functionality); the package has the functions to do the pairwise post-hoc testing, and they have an example for you to follow at the end of the documentation.

ADD COMMENTlink written 5.3 years ago by Steven Lakin1.5k
Please log in to add an answer.


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