Hi, 
I have some RNA seq data with 3 genotypes and two treatments. I'd like to do pairwise comparisons between every group. Is there an easy way to do this in ballgown, or do I need to subset my dataset into each combination of genotype/treatment, and run the ballgown pipeline on them pair by pair? What I really want as an output is a table of log2 fold change values for each pairwise comparison. 

The command where I think this may be relevent is pasted below:

> results_transcripts = stattest(bg_chrX_filt, feature="transcript", covariate="genotype_salt_mM", getFC=TRUE, meas="FPKM")
Warning message:
In stattest(bg_chrX_filt, feature = "transcript", covariate = "genotype_salt_mM",  :
  fold changes only available for 2-group comparisons

Thanks,

Rachael DeTar

Hi Rachael, ballgown doesn't do pairwise comparisons by default, but there are a few ways to force it to do it for you. The strategy you mention (running the pipeline on them pair by pair) is one option (you can use the "subset" function to make 6 ballgown objects, one with each genotype x pair treatment, then run "stattest" with each of the 6 objects). This is probably the easiest way to get what you want -- you'll have a two-group comparison and then you can pull the fold change values directly from the stattest function. (There are a couple other strategies I can think of as well, but they aren't as straightforward as the above and if you do want the 6 fold change values, there's some annoying data wrangling you'll have to do yourself. Ballgown doesn't automatically emit coefficients for more complex designs like yours because it's ~impossible to guess which values are important to the scientist in complex designs. But if this doesn't work, feel free to reply and folks can help brainstorm better solutions).