This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How to use edgeR to compare multiple conditions to reference condtiion

I have an RNA-seq experiments with a control, and 4 treatments, for each I have 3 replicates. I would like to find genes that are differentially expressed in any of these conditions compared to the control samples.

Meaning if gene A is up/down-regulated in one or more of the treatments compared to the control, I would like to know about it.

lcd = read.delim('mydata', row.names=1, as.is=T)
group = factor(c(rep("ctrl", 3), rep("treat1",3), rep("treat2",3), rep("treat3",3),rep("treat4",3)))
cds  <- DGEList(lcd, group = group)

What is the alternative for using exactTest multiple times and intersecting the results, given that I would not like to average the expression of the treatments as the different treatment have different consequences?

Thanks!

differential-expression edger rna-seq
• 17,745 views • link
updated by Ram Ram • written by xmix • xmix

1 answer 

design = model.matrix(~group)
f <- glmFit(cds, design)
lrt <- glmLRT(f, coef=2:4)
topTags(lrt)

In other words, you'd test the full model against one lacking group information.

Edit: BTW, you might find this older thread from the bioconductor list useful. It's often best to search that list when a question arises, since odds are good someone has asked a similar question there (or here) before.

• 0 views • link
written by Devon Ryan Devon Ryan

I think design = ~group should be

design = model.matrix(object = ~group)
• 0 views • link
updated by Ram Ram • written by xmix • xmix

You are correct! Good catch. I've updated the answer to correct that mistake.

• 0 views • link
written by Devon Ryan Devon Ryan

I know this is old, but why isn't it coef=2:5??

• 1 views • link
written by goodez goodez

those coefficients denote comparison of ctrl (1) with : treat (2) , treat2 (3), treat3 (4), treat4 (5). goodez

• 0 views • link
written by cpad0112 cpad0112

Yes so 4 treatments, why is there only 3 coefficients (2:4)?

• 0 views • link
written by goodez goodez

Hi, Devon.

I did not understant one thing. If I have 4 groups: 1,1,1,2,2,2,3,3,3,4,4,4, and I want to compare all of them against each other, how can I do that?

For exemplo: 

x <- read.delim("table.csv",row.names="genes")
group <- factor(c(1,1,1,2,2,2))
y <- DGEList(counts=x,group=group)
y <- calcNormFactors(y)
design <- model.matrix(~group)
y <- estimateDisp(y,design)
et <- exactTest(y)

In this case, I am comparing groups 1 and 2. Them I would do that for the other comparisons.

Or, instead of that should I do something like this:

x <- read.delim("table.csv",row.names="genes")
group <- factor(c(1,1,1,2,2,2,3,3,3,4,4,4))
y <- DGEList(counts=x,group=group)
y <- calcNormFactors(y)
design <- model.matrix(~group)
y <- estimateDisp(y,design)
fit <- glmQLFit(y,design)
qlf_quasi <- glmQLFTest(fit,coef=2:4)
topTags(qlf_quasi)

The problem is that I don't know exactly how to use coef for this kind of comparisons. coef=2:4 will compare all groups against the control group, but I do not have a control group, I want to compare everything.

Other problem is that when I use the code like this the output is a single table with multiple comparisons, but just one p-value and one FDR. I need separeted p-values for each comparison (like pos-tests for ANOVA).

It would be great if you can help me.

Thanks

• 0 views • link
written by silas008 silas008

that is also my question , i have the same situation but i couldn't reach the proper answer.

• 0 views • link
written by ali.hakimzadeh • ali.hakimzadeh

Log in to answer this question.

More posts like this