DESeq2 question regarding contrasts

disulfidebond · 10-15-2015, 05:18 PM

Hello,

I have a question regarding DESeq2 and contrasts, and I may simply be over-thinking the problem. I have a dataset with 3 groups, and two treatments

Sample Condition Group
sample1 Control Bact1
sample2 Control Bact1
sample3 Control Bact1
sample4 Treated Bact1
sample5 Treated Bact1
sample6 Treated Bact1
sample7 Control Bact2
sample8 Control Bact2
sample9 Control Bact2
sample10 Treated Bact2
sample11 Treated Bact2
sample12 Treated Bact2
sample13 Control Bact3
sample14 Control Bact3
sample15 Control Bact3
sample16 Treated Bact3
sample17 Treated Bact3
sample18 Treated Bact3

The main goal in the experiment is to determine if strain Bact1 has differentially expressed genes after treatment compared to Bact2. Previously, I did an analysis of Bact1 vs Bact2, (using the table above, minus Bact3) and was set up as:

dds <- DESeqDataSetFromHTSeqCount(sampleTable = sampleTable, directory = directory, design = ~ group + cond + group:cond)

res<- results(name="groupBact2.condTreated")

# shows interaction term, or how the treatment condition is different across Bact2 v. Bact1

res_Bact2 <- results(contrast=list(c("condition_Treated_vs_Control", "groupBact2.condTreated")))

# shows treatment condition in Bact2

# alternatively I can use dds$group <- relevel(dds$group, ref="Bact2") to look at Bact1

Now, there is a third strain in the analysis, and here is where I'm stuck. I tried creating the table above, using the design

design = ~ group + cond + group:cond

# called DESeq2() with betaprior=FALSE

and it ran successfully without errors, but the resultsNames() table was

[1] "Intercept" "group_BACT2_vs_BACT1" "group_BACT3_vs_BACT1"
[4] "cond_Trt_vs_Ctrl" "groupBACT2.condTrt" "groupBACT3.condTrt"

which doesn't seem to be what I want, because checking the attributes with

attr(dds, "ModelMatrix")

shows that the variables "group_..." doesn't include the treatment variables, and "groupBACT2.condTrt" is only the treated groups (I think?). Contrasting "groupBACT2.condTrt" and "group_BACT2_vs_BACT1" shows an output, but I'm not sure if it's doing what I'd want, i.e. contrasting the difference between Bact3 v. Bact2 in the treatment condition.

If I combine the factors into groups with

dds$grp <- factor(paste0(dds$group, dds$cond))

design(dds) <- ~ grp

running DESeq() and then resultsNames() gives:

[1] "Intercept" "grpall_Bact1Trt_vs_Bact1Ctrl"
[3] "grpall_Bact2Ctrl_vs_Bact2Ctrl" "grpall_Bact2Trt_vs_CBact1Ctrl"
[5] "grpall_Bact3Ctrl_vs_Bact1Ctrl" "grpall_Bact3Trt_vs_Bact1Ctrl"

Which is still not what I'm looking for. Am I on the wrong track entirely, or should I simply revelel() the group factors to assign the Treatment as the reference, or something else?

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
DESeq2 Contrasts	cacti	Bioinformatics	4	08-04-2015 11:20 AM
DESeq2 - using contrasts vs. multiple runs for DE analysis	rutgers2015	Bioinformatics	2	07-14-2015 10:34 AM
DESeq2 model design + contrasts	gamb	Bioinformatics	39	06-25-2015 02:06 PM
Contrasts with DESeq2	mistrm	RNA Sequencing	3	07-24-2014 09:03 AM
Differential expression analysis on RPKMs - Contrasts and contrasts of contrasts	giorgifm	Bioinformatics	6	08-16-2013 10:50 AM

10-15-2015, 05:18 PM	#1
disulfidebond Junior Member Location: United States Join Date: Oct 2015 Posts: 3	DESeq2 question regarding contrasts Hello, I have a question regarding DESeq2 and contrasts, and I may simply be over-thinking the problem. I have a dataset with 3 groups, and two treatments Sample Condition Group sample1 Control Bact1 sample2 Control Bact1 sample3 Control Bact1 sample4 Treated Bact1 sample5 Treated Bact1 sample6 Treated Bact1 sample7 Control Bact2 sample8 Control Bact2 sample9 Control Bact2 sample10 Treated Bact2 sample11 Treated Bact2 sample12 Treated Bact2 sample13 Control Bact3 sample14 Control Bact3 sample15 Control Bact3 sample16 Treated Bact3 sample17 Treated Bact3 sample18 Treated Bact3 The main goal in the experiment is to determine if strain Bact1 has differentially expressed genes after treatment compared to Bact2. Previously, I did an analysis of Bact1 vs Bact2, (using the table above, minus Bact3) and was set up as: dds <- DESeqDataSetFromHTSeqCount(sampleTable = sampleTable, directory = directory, design = ~ group + cond + group:cond) res<- results(name="groupBact2.condTreated") # shows interaction term, or how the treatment condition is different across Bact2 v. Bact1 res_Bact2 <- results(contrast=list(c("condition_Treated_vs_Control", "groupBact2.condTreated"))) # shows treatment condition in Bact2 # alternatively I can use dds$group <- relevel(dds$group, ref="Bact2") to look at Bact1 Now, there is a third strain in the analysis, and here is where I'm stuck. I tried creating the table above, using the design design = ~ group + cond + group:cond # called DESeq2() with betaprior=FALSE and it ran successfully without errors, but the resultsNames() table was [1] "Intercept" "group_BACT2_vs_BACT1" "group_BACT3_vs_BACT1" [4] "cond_Trt_vs_Ctrl" "groupBACT2.condTrt" "groupBACT3.condTrt" which doesn't seem to be what I want, because checking the attributes with attr(dds, "ModelMatrix") shows that the variables "group_..." doesn't include the treatment variables, and "groupBACT2.condTrt" is only the treated groups (I think?). Contrasting "groupBACT2.condTrt" and "group_BACT2_vs_BACT1" shows an output, but I'm not sure if it's doing what I'd want, i.e. contrasting the difference between Bact3 v. Bact2 in the treatment condition. If I combine the factors into groups with dds$grp <- factor(paste0(dds$group, dds$cond)) design(dds) <- ~ grp running DESeq() and then resultsNames() gives: [1] "Intercept" "grpall_Bact1Trt_vs_Bact1Ctrl" [3] "grpall_Bact2Ctrl_vs_Bact2Ctrl" "grpall_Bact2Trt_vs_CBact1Ctrl" [5] "grpall_Bact3Ctrl_vs_Bact1Ctrl" "grpall_Bact3Trt_vs_Bact1Ctrl" Which is still not what I'm looking for. Am I on the wrong track entirely, or should I simply revelel() the group factors to assign the Treatment as the reference, or something else?