Pathway test for reversal of effect (three-condition comparison)

Question

[glocke01]

I have an experiment with two treatments, A and B, with a 2x2 design (four conditions: ∅, A, B, and AB, with ∅ being control), and I have RNA-Seq for each condition (9 biological replicates). I want to do some sort of pathway/gene set analysis to characterize the gene expression changes that separate A from ∅ but then change in the opposite direction between AB and A. In other words, what can I learn about the effects of A monotherapy that are reversed in the dual therapy AB?

The most obvious thing is to look at the intersection of DEGs in the A - ∅ and A - AB and then do hypergeometric/overrepresentation tests against my gene sets. Along the same lines, but a bit more "rigorous" and no harder, I can do GSEA or CAMERA or such-like with both contrasts and then find the intersection of differentially expressed pathways.

However, it would be nice to do this analysis "in one step" using some sort of statistical framework that allows for three conditions in one test. Does any such test exist? I mean, my suspicion is that the answer is no and that linear models in general are not designed to do this kind of thing.

[glocke01]

It's relatively simple to describe what I'm looking for as "effects of A reversed by B", but I am aware that that language is substantively imprecise - literal reversal would be treating with A then B. The biological question of interest would be more accurately described as "effects of A - ∅ absent/attenuated in AB - ∅".

Answer 1

[alserg]

I don't know how to do it properly, but from a practical perspective you can rank genes for GSEA based on comparison (AB or zero) vs A. A signal to noise metric worked fine for me. There top genes would have similar expression levels in AB and zero, while up-regulated in A.

[glocke01]

Nice. Good thinking: rank based on some relevant comparison, send ranks to GSEA or suchlike. I guess the trouble is defining a single number for ranking purposes - by signal to noise, I take it you're thinking something like (mean in A minus mean in AB and ∅)/stdev(all samples). Something like averaging the t values (as computed by limma-voom, or stat from DESeq2) from the two comparisons would have a similar interpretation I think. (I hope - at least that's what I think those statistics are supposed to do...)

[alserg]

I got ranks = (m2 - m1) / (s1 + s2) where m2 and m1 are means, and s1 and s2 are standard deviations for the corresponding samples. This better forces top genes to have small variance between AB and zero samples.

Btw, I remembered it was used in a published paper, so it may be useful to cite: https://www.nature.com/articles/s41586-018-0052-z/figures/9

[glocke01]

I think a geometric/harmonic average might be appropriate (as opposed to the arithmetic mean implicit in my SNR statistic) so that small values in either comparison would strongly penalize the final statistic...