Hi, all. I've still been working at reconciling the differential results from this pipeline with manual DESeq2 calls. I am using the same count matrix as input and used the same model reported in the MultiQC file for the design in my DESeq2 experiment. I am doing the same contrast on the same data with the same design. However, the results I see from manually running the analysis are completely different than what the pipeline outputs.

Hi Megan. The script the pipeline runs will appear in the process directory (.command.sh). If I were you I would just compare that with whatever you are running. The pipeline also does some filtering, which may play a part. There's been quite a bit of work by other folks here to overhaul the way the modelling is handled, so things will change in the next release to allow a lot more customisation of the formula.

Is there an ETA for the next release?

not as yet

Thanks for the notes, Jonathan! I was able to work on this a bit more and am now seeing more congruence between the genes I'm calling and the ones the pipeline is calling. The gene IDs and the expression trends (i.e. treatment > control) are largely conserved but the fold change values calculated by the pipeline don't match mine. There are also a few crazy numbers reported for fold change (seen in the attached image). I'm not sure what could be causing those issues. Any ideas?

There will be a big difference between the parametric fold changes deseq2 pops out, and any empirical fold changes calculated using the ratio between groups. I forget which the pipeline outputs for that table, but suspect you're using the other!

Ahhh, okay! Thank you! That helps a lot! 🙂 I was worried when I saw such huge numbers there

Hello @Megan Justice, just wondering if you are still not managing to produce similar numbers, and if you are still getting those very large log fold changes?

Thanks for checking in! In short, I was able to produce similar results. After looking more closely at the input matrices from the two pipelines, I realized there was an error with the column order which I was able to resolve. I then reran the analysis and got much more similar numbers of differential genes. Attached is a quick graphic I made to compare the outcomes for 4 different comparisons I made from the two matrices. The results from the matrix I generated manually are labeled STAR. The core pipeline matrix results are labeled as Salmon. For comparisons 3 and 4, nearly all the DEGs I called were also called in the core pipeline. There's some disparity with comparisons 1 and 2, but the effect from those conditions is fairly small so I was not overly concerned with a lack of overlap there. As for the large fold change values shown in the image I shared previously, if you apply a log2 to them you get roughly the log2 fold change reported by my own calculation. Overall this data isn't rich with differential genes in either analysis. However, the results from the two analyses are generally similar. Thank you guys for your support in trying to figure out what was going on!