[BioC] Interaction categorical/continuous variable DESeq2

[Michael Love]

hi Hugo,

On Mon, Sep 15, 2014 at 9:41 AM, Hugo Varet <hugo.varet at pasteur.fr> wrote:

> Dear list, dear Mike Love,
>
> I am using DESeq2 to model counts from an unusual type of experiment and I
> have a question about the strategy I employed. The experiment consisted in
> sequencing 33 samples for which we have the following information:
>  - group (16 samples from group A and 17 from group B)
>  - a continuous variable X almost uniform (variable of interest)
>
> I have to add the group to the design formula because I know it has a
> strong effect on the counts. Then, as my goal is to detect genes which vary
> with the continous variable X in the same way within both groups A and B, I
> want to exclude genes for which there is an interaction between group and
> X. The design is thus ~ group + X + group:X and I used the following lines
> to test the interaction:
>
> dds <- DESeqDataSetFromMatrix(countData=counts, colData=target, design =
> ~ group + X + group:X)
> dds <- estimateSizeFactors(dds)
> dds <- estimateDispersions(dds)
> dds <- nbinomWaldTest(dds)
> res <- results(dds, name="groupB.X")
> sum(res$padj<=0.05, na.rm=TRUE)
> hist(res$padj)
>
> As I found no significant interaction (the minimum adjusted p-value is
> about 0.6), I decided to remove the interaction term from the design and to
> use ~ group + X. I can then test for the coefficients of X.
>
> If I do not detect any significant interaction, I think it is due to a
> lack of power. So, can I use the additive model ~ group + X even if it will
> not be correct for genes which actually have an interaction?
>

​I'd recommend examining the size of the estimated interaction terms in the
model​ including the interaction:

hist(res$log2FoldChange)

if these are all small, then you might prefer the simpler model, ~ group +
X.

​Just a reminder about using a continuous covariate with log link GLM: you
are looking for genes which closely follow the rule counts = 2^(a + bX).
Such a gene would follow the rule: if going from X=0 to X=0.5 produces a
fold change of 1.5, going from X=0.5 to X=1 also produces a fold change of
1.5, and X=0 to X=1 produces a fold change of 1.5^2.

​Mike​



>
> Many thanks in advance,
>
> Hugo
>
> PS: I am using R 3.1.1 and DESeq2 1.4.5
>
> _______________________________________________
> Bioconductor mailing list
> Bioconductor at r-project.org
> https://stat.ethz.ch/mailman/listinfo/bioconductor
> Search the archives: http://news.gmane.org/gmane.
> science.biology.informatics.conductor
>

	[[alternative HTML version deleted]]