I'd like to know how the log2 fold change is calculated between target and comparison population in DEXSeq. Going over the estimateExonFoldChanges function in an older version (0.12.1) of the package, I realize the interaction coefficient is taken from the model: count ~ condition * exon and fold change is calculated by applying a variance stabilizing transformation and then transformed to a log2 scale:
alleffects <- do.call(rbind, alleffects) alleffects <- vst(exp(alleffects), object) alleffects <- log2(alleffects/alleffects[, denoCol]). ###foldChange <- effects[,"target"] - effects[,"comparison"]
alleffects data frame looks like the following:
Gene target comparison Ptma 5.36425504487572 5.10532234811512 Ptma 4.43604234783272 4.7521435893567 Ptma 4.30355887270297 4.72294913039353 Tmpo 2.12202872975088 1.08346386248873 Msn 1.86941999824138 2.34083780006062
However, in a newer version that I am currently using (1.28.3), looks like the target and comparison values being used for fold change calculation are not vst transformed. Also, being divided by log(2).
alleffects <- rbind(alleffectsBM, alleffectsSM) alleffectsVst <- vst(exp(alleffects), object) alleffects <- alleffects/log(2) alleffects <- alleffects - alleffects[, denoCol] ###foldChange
The calculation in the older version makes more sense for a log2 fold change than in version 1.28.3.
An explanation for this discrepancy in fold change calculation would be appreciated.
Thanks in advance!