DeSeq design

Question

[chimerajit]

I have a simple query

when we do DSeq using two condition e.g. control and treated like this

"

Samples   Condition\
Sample1     control \
Smaple2    control\
Sample3  treated \
Sample4  treated\

And set up design parameter design =~Condition then is it analysing control vs treated or treated vs control?

Answer 1

[ATpoint]

That depends. By default control would be the reference because factors are ordered alphabetically, so it would be treated vs control, meaning that positive fold changes would indicate overexpression in treated. But you can explicitely control that. If you do res <- results(dds, contrast=c("Condition", "control", "treated")) then it would be vice versa.

> library(DESeq2)
> 
> set.seed(1)
> dds <- DESeq2::makeExampleDESeqDataSet(m=4)
> dds <- DESeq(dds)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
> 
> #/ here A is the reference by default
> dds$condition
[1] A A B B
Levels: A B
> 
> #/ that will give condition_B_vs_A by default
> resultsNames(dds)
[1] "Intercept"        "condition_B_vs_A"
> res1 <- results(dds, name="condition_B_vs_A")
> 
> #/ lets take a random gene
> gene <- rownames(res1[160,])
> 
> #/ positive fold changes means higher in B, check counts for proof
> res1[gene,]
log2 fold change (MLE): condition B vs A 
Wald test p-value: condition B vs A 
DataFrame with 1 row and 6 columns
         baseMean log2FoldChange     lfcSE      stat    pvalue      padj
        <numeric>      <numeric> <numeric> <numeric> <numeric> <numeric>
gene160   241.127        1.18405  0.412437   2.87086 0.0040936  0.999593
> 
> data.frame(condition=dds$condition,
+            counts=counts(dds, normalized=TRUE)[gene,])
        condition   counts
sample1         A 126.2945
sample2         A 168.5249
sample3         B 352.3384
sample4         B 317.3500
> 
> #/ we can turn it around:
> res2 <- results(dds, contrast=c("condition", "A", "B"))
> 
> #/ now it has a negative fold change, stats are the same
> res2[gene,]
log2 fold change (MLE): condition A vs B 
Wald test p-value: condition A vs B 
DataFrame with 1 row and 6 columns
         baseMean log2FoldChange     lfcSE      stat    pvalue      padj
        <numeric>      <numeric> <numeric> <numeric> <numeric> <numeric>
gene160   241.127       -1.18405  0.412437  -2.87086 0.0040936  0.999593

[chimerajit]

Thank for your explanation. This is a great help.