I have a ChIP-seq dataset with two conditions. For each condition I have two replicates of both input and the acutal ChIP experiment: 8 samples in total.
In R the design matrix would look like this:
myDesing <- data.frame(sample=1:8, data=c(rep('Input',4), rep('Chip',4)), condition=c('Wt','Kd'))
myDesing$condition <- factor(myDesing$condition, levels=c('Wt','Kd'))
myDesing$data <- factor(myDesing$data, levels=c('Input','Chip'))
> model.matrix(data=myDesing, ~data + condition)
(Intercept) dataChip conditionKd
1 1 0 0
2 1 0 1
3 1 0 0
4 1 0 1
5 1 1 0
6 1 1 1
7 1 1 0
8 1 1 1
I'm interested in how to handle the input data in regards to doing a differential binding analysis with edgeR, DESeq(2) or Voom-Limma.
Is it possible to incorporating the input data into the linear model?
My idea is to make a model that includes both the input and chip libraries and then test the interaction between chip and condition. More specifically for the R code example above this would correspond to testing the interaction dataChip::conditionKd. Does this make sense or does it violate some of the assumptions?
Does this not correspond to measuring the fold change of Input vs WT ( FCiwt ) and the fold change of Input vs KD (FCikd) and then comparing whether there is a difference between FCiwt and FCikd?
I have read this DESeq2 for ChIP-seq differential peaks which discuss whether one can subtract the input reads from before the DE analysis but it is a bit unsatisfactory as the discussion did not result in a good way to incorporate the input data.
In advance thanks for your help.