Question: DESeq2 3-fatcor design and different interaction terms
gravatar for c.evang
4.6 years ago by
c.evang10 wrote:


I am using DESeq2 to analyze differential expressed genes (DEGs) from counts generated by Illumina sequencing.

I have a 3-factor experiment:

  • genotype (with 3 levels): A, B, and C, each with 6 biological replicates (3 for each treatment)
  • treatment (with 2 level): mDr (the treatment) and WW (the control)
  • timepoint (with 3 levels): T1, T2 and T3.

In order to analyze genes that are influenced by an individual component and those influenced by two or all variables, I'm interested in the following results/questions:

  1. How many DEGs are influenced by an individual component (genotype, treatment, timepoint)?

  2. How many DEGs are influenced by the combination of two or all components (genotype:treatment + genotype:timepoint + treatment:timepoint + genotype:treatment:timepoint)?.

So, I tried as follows:


# Load row counts
countData<-read.delim("all_counts_matrix_T1_T2_T3", sep="\t", header=TRUE, row.names = 1)
col_ordering = c(1,3,5,2,4,6,7,9,11,8,10,12,13,15,17,14,16,18, 19,21,23,20,22,24,25,27,29,26,28,30,31,33,35,32,34,36,37,39,41,38,40,42,43,45,47,44,46,48,49,51,53,50,52,54)
countData_ord = countData[,col_ordering]

# Create countData as matrix
countData_m <- as.matrix(countData_ord)
# Interger numbers in matrix
storage.mode(countData_m) = "integer"

#Create sample information
colData<-read.delim("colData.txt", sep="\t", header=TRUE, row.names = 1)


rownames(colData) == colnames(countData_m)

# Count matrix and sample information input and design formula
dds<-DESeqDataSetFromMatrix(countData = countData_m, colData = colData, design = ~ genotype + treatment + timepoint + genotype:treatment + genotype:timepoint + treatment:timepoint + genotype:treatment:timepoint)

dds$treatment<- relevel(dds$treatment, "WW")


## Filtering to remove rows with 0 reads or only a single count across all samples
dds <- dds[ rowSums(counts(dds)) > 1, ]

# Run DESeq2
dds = DESeq(dds)
# Calling results
res <- results(dds)
mcols(res, use.names=TRUE)

# Print the coefficients

write.table(res, file='genotype_treatment_timepoint/res_04082016.txt', col.names=T, row.names=T, quote=F, sep="\t")

In resultsNames(dds) I get:

[1] "Intercept" "genotype_B_vs_A" "genotype_C_vs_A" "treatment_mDr_vs_WW"
[5] "timepoint_T2_vs_T1" "timepoint_T3_vs_T1" "genotypeB.treatmentmDr" "genotypeC.treatmentmDr"
[9] "genotypeB.timepointT2" "genotypeC.timepointT2" "genotypeB.timepointT3" "genotypeC.timepointT3"
[13] "treatmentmDr.timepointT2" "treatmentmDr.timepointT3" "genotypeB.treatmentmDr.timepointT2" "genotypeC.treatmentmDr.timepointT2" [17] "genotypeB.treatmentmDr.timepointT3" "genotypeC.treatmentmDr.timepointT3"

First of all, is it correct?

And then, how can I extract the specific results I want (question 1 and 2)? I've tried to read up on this, but haven't been able to find how to code this.

Any guidance is much appreciated.

Thanks in advance


ADD COMMENTlink modified 4.6 years ago by GenoMax96k • written 4.6 years ago by c.evang10

Hi, did you get the results what you were looking for? Would be great if you could please share how you have achieved your results. I'm currently using a similar experiment - thanks

ADD REPLYlink written 3.8 years ago by bioinfo1730
gravatar for Carlo Yague
4.6 years ago by
Carlo Yague5.7k
Carlo Yague5.7k wrote:


Look up the help page for ?results(). There are many exemple at the bottom of the page. To answer your two questions you will have to call results() several times (to test separately the impact of time points, treatments and genotypes) then combine the DEGs.

You should also read sections 1.6, 3.2, 3.3 and 3.4 from the DESeq2 manual.

Also, note that when you do the following (with your design), you only test for one contrast. Check it out using head().

res <- results(dds)
ADD COMMENTlink modified 4.6 years ago • written 4.6 years ago by Carlo Yague5.7k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 2376 users visited in the last hour