Frozen variance stabilizing transformation for count data

frozen_vst.R

mat <- matrix(rnbinom(2e5, mu=100, size=1/.01), ncol=100)

library(DESeq2)

d <- DESeqDataSetFromMatrix(mat, DataFrame(x=rep(1,100)), ~1)
# library size correction, centered log ratio to reference sample
d <- estimateSizeFactors(d)
# variance
d <- estimateDispersionsGeneEst(d)
# trend
d <- estimateDispersionsFit(d)
# approximate variance stabilizing transformation
vsd <- varianceStabilizingTransformation(d, blind=FALSE)
# geometric mean of training data for later
fixedGeoMeans <- exp(rowMeans(log(counts(d))))
# get the data
transformedData <- assay(vsd)

newmat <- matrix(rnbinom(2000, mu=100, size=1/.01), ncol=1)
d2 <- DESeqDataSetFromMatrix(newmat, DataFrame(x=rep(1,ncol(newmat))), ~1)
# re-use the reference sample
d2 <- estimateSizeFactors(d2, geoMeans=fixedGeoMeans)
# apply the trend
dispersionFunction(d2) <- dispersionFunction(d)
# compute the VST again using fixed parameters
vsd2 <- varianceStabilizingTransformation(d2, blind=FALSE)
# get the data
transformedData2 <- assay(vsd2)