Hi all,

I'd like to use DESeq2 for treatment-specific differences over time.

Here is my colData with a common time 0 for Mock, treatment1 and treatment 2:

> as.data.frame(colData(cds))
                  times    treatments
Control_a            00        Mock
Control_b            00        Mock
Control_c            00        Mock
Mock_03hr_a           03        Mock
Mock_03hr_b           03        Mock
Mock_03hr_c           03        Mock
Mock_06hr_a           06        Mock
Mock_06hr_b           06        Mock
Mock_06hr_c           06        Mock
Mock_12hr_a           12        Mock
Mock_12hr_b           12        Mock
Mock_12hr_c           12        Mock
treatment1_03hr_a      03      treatment1
treatment1_03hr_b      03      treatment1
treatment1_03hr_c      03      treatment1
treatment1_06hr_a      06      treatment1
treatment1_06hr_b      06      treatment1
treatment1_06hr_c      06      treatment1
treatment1_12hr_a      12      treatment1
treatment1_12hr_b      12      treatment1
treatment1_12hr_c      12      treatment1
treatment2_03hr_a      03      treatment2
treatment2_03hr_b      03      treatment2
treatment2_03hr_c      03      treatment2
treatment2_06hr_a      06      treatment2
treatment2_06hr_b      06      treatment2
treatment2_06hr_c      06      treatment2
treatment2_12hr_a      12      treatment2
treatment2_12hr_b      12      treatment2
treatment2_12hr_c      12      treatment2

When I tried this:

cds <- DESeqDataSetFromMatrix(countData = data, colData = colData, design = ~ treatments+times+treatments:times)

I got the error message:

Error in checkFullRank(modelMatrix) : the model matrix is not full rank, so the model cannot be fit as specified. One or more variables or interaction terms in the design formula are linear combinations of the others and must be removed.

I've been searching for the reasons and solutions like this post: With DESeq2 "Not full rank" Error with design ~ line + time + condition, and went through its comments. However, since my data is unpaired, using:

colData$times.nested = factor(c(rep(1,3),rep(rep(2:5,each=3),3)))

cds <- DESeqDataSetFromMatrix(countData = data, colData = colData, design = ~ treatments + times + treatments:times.nested + treatments:times)

or

mm = model.matrix(~ treatments + times + treatments:times, colData(cds))

cds = DESeq(cds, full=mm, betaPrior=FALSE)

gave me the same error.

How do I capture differentially expressed genes of each treatment vs. Mock over time with common time 0 incorporated?

Best regards,

Ting

hi Ting,

There was a similar question posted recently, which didn't have any treatment samples at time 0. (There was also a second level of grouping in that question, but I didn't address this in my reply to the question.) Can you take a look at my answer and then see if you have further questions:

A: deseq2 multifactor timeseries help

For building the DESeqDataSet you can use a design of ~1, and then you can use matrices for the arguments 'full' (and 'reduced' if you want to perform a likelihood ratio test).