Hi everyone,
I am trying to carry out some analysis for an experiment which has three dimensions (treatment, genotype and time) and while I have managed to extract the pairwise comparisons I am interested in, I would like to go further and look at the genes which behave differently upon my drug treatment depending on either the genotype or the time factor (for a given genotype) - I believe that this is the interaction term. However, I am not entirely sure how to do this. I have read some of the posts here, the vignette and the examples in ?results, but would like to check here.
To describe my experiment in detail: I have two genotypes (MYC On and MYC Off), which were treated with vehicle or drug (Veh and Drug), and the mRNA levels measured after 2, 6 and 18 hours of exposure.
To obtain the pairwise comparisons, I set the design using a grouping variable, and then extracted each comparison of interest:
dds <- DESeqDataSetFromMatrix(countData = countdata, colData = coldata, design = ~Genotype + Time + Treatment)
dds$group <- factor(paste0(dds$Genotype, dds$Time, dds$Treatment))
design(dds) <- ~group
#Extract each pairwise comparison, e.g.
MYC_ON_18h_TreatedvsWT <- lfcShrink(dds, contrast=c("group", "ONT18DRUG", "ONT18VEH"), type="ashr")
This worked well, but as stated, I would now like to see which which genes are differentially affected by my treatment based on the:
a) Genotype (for each timepoint, if possible).
b) Time (for a given genotype).
I am unsure as to how to write the code to pull this information out. I am under the impression that I cannot do this anymore using a grouping variable, and have to define the design formula explicitly, based on the phrasing on ?results. I am not sure how then to extract the relevant data for each timepoint and comparison. I have written something which I think might work, but would appreciate another opinion.
dds <- DESeqDataSetFromMatrix(countData = countdata, colData = coldata, design = ~Genotype + Time + Treatment + Genotype:Treatment + Time:Treatment)
#check reference levels
dds$Genotype
#Answer: Levels: ON OFF
dds$Time
#Answer: T18 T2 T6
dds$Treatment
#Answer: DMSO DRUG
#Therefore current reference level is MYC ON, vehicle, 18h.
#It would make more sense to me from an experimental flow standpoint if the current reference level is MYC ON, vehicle, 2h, so I change this below.
dds$Time = relevel(dds$Time, "T2")
dds <- DESeq(dds)
#Check the comparisons made
resultsNames(dds)
#[1] "Intercept" "Genotype_OFF_vs_ON" "Time_T18_vs_T2" "Time_T6_vs_T2"
#[5] "Treatment_DRUG_vs_DMSO" "Genotype_OFF.TreatmentDRUG" "TimeT18.TreatmentDRUG" "TimeT6.TreatmentDRUG"
#Extract results for 2h timepoint - differential effect of drug based on genotype.
results2h <- results(dds, name="Genotype_OFF.TreatmentDRUG")
#Extract results for 6h timepoint - differential effect of drug based on genotype.
#Extract results for 18h timepoint - differential effect of drug based on genotype.
#Extract results for MYC on - differential effect of drug based on time.
2v6 <- results(dds, name=" TimeT6.TreatmentDRUG") #(2h vs 6h, MYC ON)
2v18 <- results(dds, name=" TimeT18.TreatmentDRUG") #(2h vs 18h, MYC ON)
6v18 <- results(dds, contrast=list("TimeT18.TreatmentDRUG", "TimeT6.TreatmentDRUG") #(6h vs 18h, MYC ON).
#Extract results for MYC off - differential effect of drug based on time.
However, I am not sure how to extract the resutls for the 6h timepoint - differential effect of drug based on genotype, 18h timepoint - differential effect of drug based on genotype or the results for MYC off - differential effect of drug based on time. Would I have to re-define the reference level and then rerun dds? Example below:
dds$Time = relevel(dds$Time, "T6")
dds <- DESeq(dds)
#Extract results for 6h timepoint - differential effect of drug based on genotype.
results6h = results(dds, name="Genotype_OFF.TreatmentDRUG")
dds$Time = relevel(dds$Time, "T18")
dds <- DESeq(dds)
#Extract results for 18h timepoint - differential effect of drug based on genotype.
results18h = results(dds, name="Genotype_OFF.TreatmentDRUG")
Any help would be greatly appreciated.
