Interaction analysis using DESeq while adjusting for covariates.

Question

[Sumit Paliwal]

Hi All,

I have performed an RNA-Seq wherein I want to look for differential expression based on two variables: Condition (high and low) and Treatment (A and B) as well as their interaction. At the same time, I have to adjust for covariates such as 'Batch' and 'Gender'. To do this, I have specified the design as following:

design <- ~Batch + Gender + Condition + Treatment + Condition:Treatment
resultsNames(dds) 
[1] "Intercept"           "Batch_2_vs_1"       "Batch_3_vs_1"        "Gender_2_vs_1"         "Condition_high_vs_low" 
[6] "Treatment_B_vs_A"   "Conditionhigh.TreatmentB"

Can anyone suggest, if this is the right design? If so, what contrast will give me 'Batch' and 'Gender' adjusted results.

Thanks in advance.

Sumit

Answer 1

[Devon Ryan]

Your design looks reasonable and you don't need any special contrast to get the adjusted results, since you can ONLY get the adjusted results. That is, the Treatment_B_vs_A coefficiant will directly give you the results of treatment A vs. B while adjusting for batch and gender. The same holds for the other contrasts and coefficients. That's the benefit of putting batch and gender directly into the design like you wisely did, you've now made your life notably easier.

[Sumit Paliwal]

Hi Devon, Thanks for your response. I now feel more confident about my approach. Based on the above design, I am analysing my result as follows:

# The effect of Treatment on Low-Condition group 
res1 = results(dds, contrast = c("Treatment","B","A"))

# The effect of Treatment on High-Condition group 
res2 = results(dds, list(c("Treatment_B_vs_A"," Conditionhigh.TreatmentB")))

#  The difference between High-Condition and Low-Condition group without Treatment 
res3 = results(dds, contrast =c("Condition","high","low"))

# With Treatment, the difference between High-Condition and Low-Condition group 
res4 = results(dds, list(c("Condition_high_vs_low"," Conditionhigh.TreatmentB ")))

For Condition 'low' and for Treatment 'A' are my reference levels. Do you agree with this interpretation?

[swbarnes2]

Personally, I think it's more readable to do these particular contrasts by making a Condition_Treatment column, and using that.

[Devon Ryan]

We would normally interpret these as follows:

• res1: The effect of Treatment. You could just use name="Treatment_B_vs_A" instead.
• res2: Yes it's the effect of treatment on the high condition group, but I question whether this is useful. You likely want the interaction instead (name="Conditionhigh.TreatmentB")
• res3: You stated this nicely. You can also simply use name="Condiction_high_vs_low"
• res4: Like with res2, I'm not sure how relevant of a result this is.

In general with a 2x2 design you're interested in the 2 "main effects" (Treatment and Condition here) and the interaction between them (i.e., any unexpected change you get when they're combined). For most people I suspect it's simpler to use the coefficient names that resultsNames() returns.

[Sumit Paliwal]

Devon, Thanks again for your time and effort. I did not mention it earlier but I have also done (name="Conditionhigh.TreatmentB"). In the context of our study design, Yes the 'interaction' is important but the other results are also equally important. I just need to find the best possible explanation for the interaction results.

[Thomas Neff]

Hi,

I will advice you to use the makeContrasts function. You will find more details by following these links

R documentation about limma::makeContrasts

edgeR user guide check specific experimental design chapter. It is the same way with DESeq2 method

For the conditionnal tests design, it would be easier to create a new column in metadata before to perform the test. This methode is also explain in edgeR user guide.

A good exercise can be to use Likelihood Ratio test in DESeq2 method. This will force you to deepen a little into the concept of design because it is necessary to create full and reduce design in this case.

And you can check your results with edgeR, which should be similar (correlated p values)

[Sumit Paliwal]

Thanks Thomas. I guess, you mean the Likelihood Ratio test in DESeq2 and not the Log Ratio test. Will definitely have a look and try to implement your suggestions as well.

[Thomas Neff]

Yes, I was tired sorry. I will correct my message thank you !