I have RNAseq data from 5 patients before and after treatment, and 5 controls before and after treatment. I would like to identify genes that respond differently to treatment in patients compared to controls. Because we have data from the same individuals before and after treatment, we would like the model to include a term for subjects, ie. design=~Disease+Disease:Subject+Disease:Treatment.
While this comparison is possible using EdgeR, I cannot get it to work using DESeq2. I keep receiving the message "Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) : contrasts can be applied only to factors with 2 or more levels". Is there any way to make this comparison using DESeq2?
My script and the error I receive is below:
> coldata = with(samples, data.frame(LibraryName = I(LibraryName),Disease = Disease,Treatment = Treatment,Subject = Subject))
> dds<-DESeqDataSetFromMatrix(countData=countData, colData=coldata, design=~Disease+Disease:Subject+Disease:Treatment)
> dds$Disease<-relevel(dds$Disease, "HAP")
> dds$Treatment<-relevel(dds$Treatment, "Control")
estimating size factors
gene-wise dispersion estimates
final dispersion estimates
fitting model and testing
Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
contrasts can be applied only to factors with 2 or more levels
Thank you in advance for your suggestions!