lfcshrink in DESeq2 with grouping variable
1
0
Entering edit mode
kendralec ▴ 10
@kendralec-14084
Last seen 6.3 years ago
Stanford University

Hello,

I am looking for advice about how to correctly use the "coef" argument to lfcshrink when I am using DESeq2 with a grouping variable. I am analyzing an RNAseq dataset with multiple factors and variables using DESeq2:

  • SurgeryGroup (injury versus sham)
  • TimePoint (1 versus 2)
  • Sex (M versus F)

My exploratory analysis of the data indicated that there is (surprisingly) very little effect of Sex so I control for it in my model and focus on the other variables. I am interested in asking the following questions:

  1. What genes are differentially expressed in injury versus sham at timepoint 1?

  2. What genes are differentially expressed in injury versus sham at timepoint 2?

To this end, I have created a grouping variable with levels corresponding to the 4 possible combinations of SurgeryGroup and TimePoint. I can then call results and shrink the LFCs by specifying my contrasts of interest.

> dds$group <- factor(paste0(dds$SurgeryGroup, dds$TimePoint))
> design(dds) <- ~ Sex + group
> dds <- DESeq(dds)
> resultsNames(dds)
[1] "Intercept"                "Sex_M_vs_F"               "group_sham1_vs_sham2"
[4] "group_injury2_vs_sham2" "group_injury1_vs_sham2"

> resultsTimePoint1 <- results(dds, contrast=c("group", "injury1", "sham1"), alpha=0.05)
> resultsTimePoint1 <- lfcShrink(dds, contrast=c("group", "injury1", "sham1"), res=resultsTimePoint1)

However, I'm interested in using "apeglm" to shrink the LFCs, which requires the "coef" argument to be used instead of "contrast" for lfcshrink. Is it possible to do this with a grouping variable as I have done? Perhaps I am just confused about how to specify the correct contrasts I am interested in by its coefficient, so if that's the case, I'd appreciate advice about how to do this.

Thanks for your consideration!

deseq2 lfcshrink • 1.3k views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 6 hours ago
United States

Here is some example code.

You'd just add sex to the beginning of the design to adapt to your code, ~sex + time + time:condition.

> dds <- makeExampleDESeqDataSet()
> dds$condition
 [1] A A A A A A B B B B B B
Levels: A B
> dds$time <- factor(rep(c(1,2,1,2),each=3))

> design(dds) <- ~time + time:condition

> dds <- DESeq(dds, quiet=TRUE)
> resultsNames(dds)
[1] "Intercept"        "time_2_vs_1"      "time1.conditionB" "time2.conditionB"

> res1 <- lfcShrink(dds, coef=3, type="apeglm")
> res2 <- lfcShrink(dds, coef=4, type="apeglm")

These last two could equivalently have been:

> res1 <- lfcShrink(dds, coef="time1.conditionB", type="apeglm")
> res2 <- lfcShrink(dds, coef="time2.conditionB", type="apeglm")

These two are the condition B vs A effect for time 1 and for time 2, respectively.

 

ADD COMMENT
0
Entering edit mode

Hi Michael, Thanks for your quick response! This works perfectly. I was having trouble wrapping my brain around the examples given in the "Group-specific condition effects" section of the vignette, but it makes sense now. Much appreciated!

 

ADD REPLY

Login before adding your answer.

Traffic: 620 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6