DESeq2 design for differential expression with 2 timepoints and 2 controls

Question

[Cristina]

Hi everyone! I’m trying to analyze RNA-seq data using DESeq2 to perform differential expression analysis, but I’m not sure if I’m using the correct design for this type of experiment or if I’m extracting the results correctly. I searched for similar experimental designs but couldn’t find exactly what I needed.

To give you some context: We supplemented a bacterial strain with a fatty acid (Fed condition). As a control, we used cultures supplemented with DMSO (the solvent used to dissolve the fatty acid) to ensure that any changes are due to the fatty acid itself and not the solvent. We harvested cells at two time points: after 30 minutes and after 6 hours of feeding. So, our experiment looks like this:

| sampleName | condition | timepoint |
| ---------- | --------- | --------- |
| DMSO1T1    | DMSO      | 30min     |
| DMSO2T1    | DMSO      | 30min     |
| DMSO3T1    | DMSO      | 30min     |
| Fed1T1     | Fed       | 30min     |
| Fed2T1     | Fed       | 30min     |
| Fed3T1     | Fed       | 30min     |
| DMSO1T2    | DMSO      | 6h        |
| DMSO2T2    | DMSO      | 6h        |
| DMSO3T2    | DMSO      | 6h        |
| Fed1T2     | Fed       | 6h        |
| Fed2T2     | Fed       | 6h        |
| Fed3T2     | Fed       | 6h        |

Since we want to study the effects of both time and feeding, I used this design:

design = ~ timepoint + condition + timepoint:condition

    dds <- DESeqDataSetFromMatrix(
  countData = cts[, sampleName],  
  colData = colData,
  design = ~ timepoint + condition + timepoint:condition
)

To answer my biological questions, I extracted results like this:

1) Do genes change over time (30min vs 6h), regardless of feeding?

res_time_30vs6 <- results(dds, name = "timepoint_6h_vs_30min")

2) How do genes change after 30 min of feeding (Fed vs DMSO at 30 min)?

res_feeding_30min <- results(dds, name = "condition_Fed_vs_DMSO")

3) How do genes change after 6 h of feeding (Fed vs DMSO at 6 h)?

res_feeding_6h <- results(dds, contrast = list(
      c("condition_Fed_vs_DMSO", "timepoint6h.conditionFed")
    ))

4) Is the effect of feeding different at 6 h compared to 30 min for Fed samples? (interaction)

res_time_fed <- results(dds, contrast = list(c("timepoint_6h_vs_30min", "timepoint6h.conditionFed")))

However, I’m still confused — I’m not sure I fully understand how the interaction works in the design, and I’m worried I didn’t specify the design or contrasts correctly. Can someone help clarify this?

Thank you so much in advance!

[jkim]

Have a look at the contrast and interaction section.

The key point to remember about designs with interaction terms is that, unlike for a design ~genotype + condition, where the condition effect represents the overall effect controlling for differences due to genotype, by adding genotype:condition, the main condition effect only represents the effect of condition for the reference level of genotype (I, or whichever level was defined by the user as the reference level). The interaction terms genotypeII.conditionB and genotypeIII.conditionB give the difference between the condition effect for a given genotype and the condition effect for the reference genotype.

Q1) is difference between 6h and 30min at the reference level of the condition since the interaction term is included in the model.

[Cristina]

Thank you so much for the clarification! I hadn’t understood that when I read about the design. With the interaction term, the main effects (like condition_Fed_vs_DMSO or timepoint_6h_vs_30min) are not global effects — they’re interpreted at the reference level of the other factor. Now i got it, thank you!

Answer 1

[swbarnes2]

Adding interactions is not correct for most of what you are doing. Time + treatment is what you want for 1. For questions like 2 and 3, make a new column of colData that is time_condition, and use that to compare subsets of samples to each other.

For 4, do use the interaction, use resultsNames to get the exact name of what you want, it should be the last option.

[Cristina]

Thanks so much for your response. I had assumed I could test main effects and interactions directly, but now I see how including an interaction changes the interpretation of the main effects. Also, your suggestion about creating a new column for specific condition-time combinations makes perfect sense. Really appreciate your help! Thank you!

[txema.heredia]

As far as I understand interaction designs, timepoint6h.conditionFed should show you if there is any synergistic effect between diet and time that is different from just the linear changes of a change in diet at t=30min and + changes in DSMO between 30min and 6h.

I always begin the analysis with a simple design model using a grouping of all condition variables (or additive if I have other covariates like the mouse's sex, or batch effects). And after that, start toying with interaction designs to see if the observed changes are just linear/additive or if there are some genes "changing more than what can be explained by linear effects".

Keep in mind that you might find many genes deemed as significant in the interaction model that actually aren't that important. For example, if both treatment and time downregulate a gene to almost 0 counts (say, from 100 in controls to 3 in diet/time = -97), the interaction design will say that such gene is upregulated because the real data/biological limitations cannot go as low as the model is predicting.

[SamGG]

Having a similar design and a new column of colData that is time_condition, could you confirm that for 4, the model with interaction leads to the same result as the following contrast: (Fed_6h - DMSO_6h) - (Fed_30min - DMSO_30min) which is easier for my mind?