Question: lfc shrinkage two factor experiment
0
5 months ago by
rbutler0
rbutler0 wrote:

Hello, I am working through the vignettes and examples for DESeq2 ?results. I had a question(s) regarding lfcshrinkage. Looking at Example 2:

 ## Example 2: two conditions, two genotypes, with an interaction term

dds <- makeExampleDESeqDataSet(n=100,m=12)
dds\$genotype <- factor(rep(rep(c("I","II"),each=3),2))

design(dds) <- ~ genotype + condition + genotype:condition
dds <- DESeq(dds)
resultsNames(dds)

# Note: design with interactions terms by default have betaPrior=FALSE

# the condition effect for genotype I (the main effect)
results(dds, contrast=c("condition","B","A"))

# the condition effect for genotype II
# this is, by definition, the main effect *plus* the interaction term
# (the extra condition effect in genotype II compared to genotype I).
results(dds, list( c("condition_B_vs_A","genotypeII.conditionB") ))

# the interaction term, answering: is the condition effect *different* across genotypes?
results(dds, name="genotypeII.conditionB")


How would you use the lfcShrink() function in those three cases (in particular with 'apeglm', which won't work with contrast in the second results() call)? Are you supposed to apply shrinkage to the interaction term? How would you specify a different alpha (0.05)?

deseq2 lfcshrink • 136 views
modified 5 months ago by Michael Love25k • written 5 months ago by rbutler0
Answer: lfc shrinkage two factor experiment
2
5 months ago by
Michael Love25k
United States
Michael Love25k wrote:

Have you seen the extended section on shrinkage methods in the vignette? apeglm only can be used for individual terms however we show some examples for how to rearrange terms to apply apeglm. Alternatively you can use ashr. You can change alpha by running results() first and supplying res=res.

ADD COMMENTlink modified 5 months ago • written 5 months ago by Michael Love25k

Thanks! I saw the refactoring and I think I understand it, but I am curious about interpreting the outcome. I can for instance relevel genotype on II and then it will become a single coef term, but then would the condition effect for this re-leveled genotype II still be the sum of the condition effect for genotype I + interaction term in the original model?

For context, we are looking at a situation where a treatment had an expected behavioral effect only in males. So ultimately the goal is to look at DE genes in the males, but separate that gene set into groups for which the condition effect is different across sex, and those for which it is not. So with a ~sex + treatment + sex:treatment model, that would be "treatment_B_vs_A" DE genes, separated by whether they are also DE in "sexF.treatmentB".

Correct? or should I be refactoring this into a single grp variable and looking at something like grp_Mtrt_vs_Muntrt to get my DE gene set and then splitting them by whether they also are in the DE gene set grp_Ftrt_vs_Mtrt?

1

This sound like the refactoring would be possible but introduce lots of placing where coding mistakes could happen. Why not just use ashr? It performs very well.