I have RNA-seq data from 2 different genotypes, and each includes a treated versus mock-control. I'd like to use the hisat2+stringtie+ballgown pipeline to identify D.E. genes between the two genotypes, in treated vs control.

I am unable to come up with the right stattest argument to code for this.

For e.g. if my phenotype file is as follows: (design does include 3x replicates)

SAMPLE(name) GENOTYPE (KO or WT) TREATMENT (COMPOUND or MOCK)

And I use stattest (bg, feature = "transcript", covariate = "GENOTYPE", adjustvars = c("TREATMENT"), getFC = TRUE, meas = "FPKM") I'm afraid I'm going to get all the transcripts that are D.E in KO vs WT, regardless of whether they have been treated with Compound or mock.

Whereas what I want is all transcripts D.E. in KO vs WT, specifically when treated with COMPOUND (not when treated with MOCK).

Thanks in advance for reading,