So recently we have got sample back with the following design
data.frame(Condition=rep(rep(c("Case","Control"), each=3),2), Treatment=rep(c("Treated", "Untreated"),each=6), Lane=c(1,2,1,2,1,2,1,2,2,2,1,1), Batch=c(1,1,4,4,1,1,4,3,2,1,3,2), RIN=c(7.7,7.7,7.6,8.1,7.8,7.7,7.5,7.8,7.9,7.4,8,7.8))
What we would like to do is to:
1. Compare effect of treatment in case
2. Compare effect of treatment in control
3. Compare and contrast the effect of treatment in case when compared to control.
From previous answers, it seems like there are two different model that we can use:
Where group is <Treatment><Condition>
Take for example, when we try to compare the difference of condition when there is no treatment
e.g. Untreated Case vs Untreated Control
For 1, we use
dds = DESeqDataSetFromMatrix(countData=countdata, colData=coldata, design = ~Batch+RIN+Lane+Condition+Condition:Treatment) dds$Condition = relevel(dds$Condition, ref="Control") dds$Treatment = relevel(dds$Treatment, ref="Untreated") dds = DESeq(dds) results(dds, contrast=c("Condition", "Control", "Case"))
For 2 we use
coldata$Group = paste0(coldata$Treatment, "-", coldata$Condition) dds = DESeqDataSetFromMatrix(countData=countdata, colData=coldata, design = ~Batch+RIN+Lane+Group) dds = DESeq(dds) results(dds, contrast=c("Group", "Untreated-Control", "Untreated-Case"))
However, I found that the two design actually give very different results.
Using the interaction term model, we will get around 438 genes with padj < 0.05 whereas none of the genes can anywhere close to significance when we use the group model.
Is there something that I did wrongly? Which one is the more appropriate method? And most importantly, why was there such a large difference?