Hi,

I still have a few questions regarding the contrast and interaction in DESeq2 after going through the forum.

1. Test on condition effect: the difference between using group variable and interactions.

There seems two ways to get the condition effect for one set basing on  the examples parts of results:

# first one:
# the condition effect in set Z
# (the interaction effect added to the main effect)
dds <- makeExampleDESeqDataSet(n=500, m = 18, betaSD = 4)
dds$set <- factor(rep(rep(c("X","Y","Z"),each=3),2))
design(dds) <- ~ set + condition + set:condition
dds <- DESeq(dds)
res <- results(dds, contrast=list(
        c("conditionB","setZ.conditionB"),
        c("conditionA","setZ.conditionA")))

# second one:
# using a grouping variable.

dds.g <- dds
dds.g$group <- factor(paste0(dds.g$set, dds.g$condition))
design(dds.g) <- ~ group
dds.g <- DESeq(dds.g)
res.g <- results(dds.g, contrast=c("group","ZB","ZA"))


cor(res[, 2], res.g[, 2])

[1] 0.9873791

The first approach seems provides larger abs(fold change) and found more DE genes.

mean(abs(res[, 2]) - abs(res.g[, 2]))

0.3428656

Number of DE genes in first versus second approach:
285 vs 257

My questions are:

(A) Are I correct in understanding that two approaches are basically doing the same thing?
(B) Which approach is more reasonable?

2.  Including betaPrior in the model.

It is mentioned  by Mike several times that if there are interaction terms in the model, betaPrior is recommended to be set FALSE.

dds.nb <- dds
dds.nb <- DESeq(dds.nb, betaPrior = FALSE)

If I am interested in genes show difference between A and B aing, then command is

res.nb <- results(dds.nb, contrast=list(c("condition_B_vs_A","setY.conditionB", "setZ.conditionB")))

(C) is it right?

3. Say I am interested in finding genes show interaction between conditions in set Z .

From the example of results,

res.i <- results(dds, contrast=list("setZ.conditionB","setZ.conditionA"),  independentFiltering = FALSE)

In the condition of betaPrior = FALSE,
 

resultsNames(dds.nb) gives
 "Intercept"        "set_Y_vs_X"       "set_Z_vs_X"       "condition_B_vs_A" "setY.conditionB"  "setZ.conditionB"

then "setZ.conditionB"  should be interactive terms, right?

res.nb.i <- results(dds.nb, name="setZ.conditionB", independentFiltering = FALSE)

However, the fold change is quite different compared to res.i.

So, (D) am I correct to calculate res.nb.i and (E) is the difference between fold change caused  by diabling of shrinkage?

Thanks in advance for your time and patience.

A,B: See my answer to a recent question:

A: DESeq2 results example

C: No, that comparison doesn't look reasonable, but I don't know what exactly you want to test. You need to be more specific.

D: yes, and E: the difference in LFC is due to shrinkage.