I am analyzing my own RNA-seq data for the first time using the Pertea Nature Protocol paper as a reference. I think ballgown is the best option for my analysis because my data is multivariate. 

My data has 3 main multivariate conditions. Group = control and experimental, Region = cortex and striatum, and treatment = control and treated. (n=3 per group for a total of 24 samples)

I can use the stattest to get the differential expression between each variable while adjusting for the other two variables. Using this one for example:

Group_genes = stattest(bulk_mg_filt, feature="gene", covariate="Group", adjustvars = c("Region","Trt"), getFC=TRUE, meas="FPKM”)

However, (if i understand correctly) this will give me the differential expression that is mostly concordant across the adjusted variables. How do I find differential expression between the group variable that maybe discordant between the adjusted variables? For example, things that are specific to one region, opposite between regions, or are only revealed after treatment and not in control treatment? Is there a way to specifically look for those discordant differences between the 2 groups? 

Thanks, Lindsay