Time course RNA seq with DESEq2

Question

[synat.keam]

Hi all,

Hope you are well. I am running time course experiment to examine differentially expressed genes from tumour tissues between two conditions, radiotherapy and sham-radiotherpapy groups at 4 time points post radiotherapy. Each time point had size matched untreated controls. I used the following code to run model with interaction.

Option 1

dds_all <- DESeqDataSetFromMatrix(countData = Raw.Count, 
                               colData = Phenotype,
                               design = ~ Treatment + Time + Treatment:Time)

set reference group

dds_all$Treatment <- relevel(dds_all$Treatment, ref ="Sham-radiotherapy")

DESeq

dds_all <- DESeq(dds_all)

Examine result

resultsNames(dds_all)

all_time_interaction <- results(dds_all, alpha =0.05)

summary(all_time_interaction)

with this model, I got around 100 differentially expressed genes (not many)

Additive model

dds_additive <- DESeqDataSetFromMatrix(countData = Raw.Count, 
                                   colData = Phenotype,
                                   design = ~ Treatment + Time)

relevel

dds_additive$Treatment <- relevel(dds_additive$Treatment, ref ="Sham-radiotherapy")

DESeq

dds_additive <- DESeq(dds_additive)

Examine result

resultsNames(dds_additive) 

all_time_additive <- results(dds_additive,alpha = 0.05)

summary(all_time_additive)

Additive model (no interaction) gave me around 3000 differentially expresse genes.

Hard to make a decision whether to retain the interaction or drop it because I did not see any p-value for interaction term? Does anyone knows where I can I look whether interaction is signficant? so I know whether I should keep it in the model or not? what is your suggestion to move forward!

Option 2

I found tutorial using the following command to analyse Time course RNA seq using DEseqDataset function. I then tried it with my dataset using following command.

Time_course_RT <- DESeqDataSet(Raw.Count, ~ Treatment + Time + Treament:Time)

Time_course_RT <- DESeq(Time_course_RT, test="LRT", reduced = ~ Treatment + Time)

My question is whether my option 1 using DESeqDataSetFromMatrix is correct to analyse Time course RNA seq compared with DEseqDataSet?

I got error once I used DEseqDataSet with the following message

*Error in DESeqDataSet(Raw.Count, ~Treatment + Time + Treament:Time) : 
  'se' must be a RangedSummarizedExperiment object*. 

How do I fixed it? how could I convert my data into RangedSummarizedExperiment object ?

I found some people experienced this issue too, but did not get any definitive answer. Thank in advance for your help.

Regards,

Synat,

Answer 1

[swbarnes2]

Running results without specifying a contrast or a ResultsName is dangerous, are you absolutely sure you know what results you are getting?

[synat.keam]

HI Swbarnes,

This is time course examining genes that are differentially expressed overtime, not at specific time point. So, in this case, I have only two groups, radiotherapy and sham-radiotherapy collecting at different time points and Time variable is now included in the model and also Time:Treatment interaction. I have specified reference group which is sham-radiotherapy so the comparision is between radiotherapy vs. sham-radiotherapy adjusting for time or time:treatment interaction. I used to specified contract previously, but setting reference level is enough for two groups. I have checked already and the model did compare between radiotherapy_vs_sham-radiotherapy. Hope i am right!

Thanks, Synat