Dear all,

although I read through many advices and posts concerning the lfc Shrink function, I am still insecure if I used it correctly. I have 4 different conditions in my RNA seq data: CTRL_NoTreatment, CTRL_Treatment, Knockdown(KD)_NoTreatment and Knockdown_Treatment and I would like to do the following comparisons:

condition_CTRL_Treatment_vs_CTRL_NoTreatment, condition_KD_NoTreatment_vs_CTRL_NoTreatment, condition_KD_Treatment_vs_KD_NoTreatment and condition_KD_Treatment_vs_CTRL_Treatment.

Is it correct to run the dds function several times to obtain all required coef I need using different reference levels?

sampleTable$condition <- relevel(sampleTable$condition, ref = "CTRL_NoTreatment")
ddsHTSeq <- DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory = read_dir, design= ~ condition)
dds <- DESeq(ddsHTSeq)
resLFC1 <- lfcShrink(dds, coef= "condition_CTRL_Treatment_vs_CTRL_NoTreatment", type = "apeglm")
resLFC2 <- lfcShrink(dds, coef= "condition_KD_NoTreatment_vs_CTRL_NoTreatment", type = "apeglm")

sampleTable$condition <- relevel(sampleTable$condition, ref = "KD_NoTreatment")             
ddsHTSeq <- DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory = read_dir, design= ~ condition)
dds <- DESeq(ddsHTSeq)
resLFC3 <- lfcShrink(dds, coef = "condition_KD_Treatment_vs_KD_NoTreatment", type = "apeglm")

sampleTable$condition <- relevel(sampleTable$condition, ref = "CTRL_Treatment")  
ddsHTSeq <- DESeqDataSetFromHTSeqCount(sampleTable=sampleTable, directory = read_dir, design= ~ condition)
dds <- DESeq(ddsHTSeq)
resLFC4 <- lfcShrink(dds, coef = "condition_KD_Treatment_vs_CTRL_Treatment", type = "apeglm")

I'm using R version and DESeq2 version 1.28.1.

Thanks a lot in advance to help with my confusion and kind regards,

Lissy

You actually only need to run nbinomWaldTest after re-leveling. So you run DESeq() once, then re-level, then run nbinomWaldTest (to form the MLE estimates) and lfcShrink.

See here:

In this case, the dispersion does not have to be estimated again, as the designs are equivalent, up to the meaning of the coefficients. Instead, one need only run nbinomWaldTest to re-estimate MLE coefficients

https://bioconductor.org/packages/release/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#extended-section-on-shrinkage-estimators