Question: DESeq2: Difference in rlog and vst values if run transform before or after running DESeq
gravatar for Zach Roe
2.8 years ago by
Zach Roe10
Zach Roe10 wrote:


I notice that the rlog and vst values for my assay differ depending if I run the DESeq function before or after running the transformation.

I thought that did this not matter as the manual states under section Data transformations and visualization -- "However, the running times are shorter and more similar with blind=FALSE and if the function DESeq has already been run, because then it is not necessary to re-estimate the dispersion values," so I assume the only difference was running time, and the transforms are performing the estimation of dispersion values if DESeq was not run prior?

Should DESeq be run before transforms?

# Case 1

dds1 <- DESeqDataSetFromMatrix(countData = countdata, colData = coldata, design = ~ condition)
rld1 <- rlog(dds1, blind=FALSE)
vsd1 <- varianceStabilizingTransformation(dds1, blind=FALSE)


# Case 2

dds2 <- DESeqDataSetFromMatrix(countData = countdata, colData = coldata, design = ~ condition)
dds2 <- DESeq(dds2)
rld2 <- rlog(dds2, blind=FALSE)
vsd2 <- varianceStabilizingTransformation(dds2, blind=FALSE)

length(which((assay(dds1) == assay(dd2)) == FALSE))
# [1] 0
length(which((assay(rld1) == assay(rld2)) == FALSE))
# [1] 583555 # the entire set
length(which((assay(vsd1) == assay(vsd2)) == FALSE))
# [1] 583555 # the entire set

I am running DESeq2_1.14.1 

R version 3.3.2 (2016-10-31)
Platform: x86_64-apple-darwin13.4.0 (64-bit)
Running under: macOS Sierra 10.12.3

ADD COMMENTlink modified 2.8 years ago by Michael Love26k • written 2.8 years ago by Zach Roe10
Answer: DESeq2: Difference in rlog and vst values if run transform before or after runni
gravatar for Michael Love
2.8 years ago by
Michael Love26k
United States
Michael Love26k wrote:

Not sure why you're not getting equality, at least within tolerance. Can you use all.equal(x,y) instead, so we can see what the mean relative difference is?

Here's what I get:

> dds <- makeExampleDESeqDataSet()
> rld <- rlog(dds, blind=FALSE)
> dds <- DESeq(dds)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
> rld2 <- rlog(dds, blind=FALSE)
> all.equal(assay(rld), assay(rld2))

[1] TRUE

> sessionInfo()
R version 3.3.2 (2016-10-31)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 16.10

 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=en_US.UTF-8       
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                  LC_ADDRESS=C              

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices datasets  utils     methods   base     

other attached packages:
 [1] DESeq2_1.14.1              SummarizedExperiment_1.4.0 Biobase_2.34.0            
 [4] GenomicRanges_1.26.1       GenomeInfoDb_1.10.1        IRanges_2.8.1             
 [7] S4Vectors_0.12.0           BiocGenerics_0.20.0        rmarkdown_1.2             
[10] magrittr_1.5               knitr_1.15.1               devtools_1.12.0           
[13] BiocInstaller_1.24.0      

loaded via a namespace (and not attached):
 [1] genefilter_1.56.0    locfit_1.5-9.1       splines_3.3.2        lattice_0.20-34     
 [5] colorspace_1.3-1     htmltools_0.3.5      survival_2.40-1      XML_3.98-1.5        
 [9] foreign_0.8-67       withr_1.0.2          DBI_0.5-1            BiocParallel_1.8.1  
[13] RColorBrewer_1.1-2   plyr_1.8.4           stringr_1.1.0        zlibbioc_1.20.0     
[17] munsell_0.4.3        gtable_0.2.0         memoise_1.0.0        evaluate_0.10       
[21] latticeExtra_0.6-28  geneplotter_1.52.0   AnnotationDbi_1.36.0 htmlTable_1.7       
[25] Rcpp_0.12.8          acepack_1.4.1        xtable_1.8-2         scales_0.4.1        
[29] backports_1.0.4      Hmisc_4.0-0          annotate_1.52.0      XVector_0.14.0      
[33] gridExtra_2.2.1      ggplot2_2.2.0        digest_0.6.10        stringi_1.1.2       
[37] grid_3.3.2           rprojroot_1.1        tools_3.3.2          bitops_1.0-6        
[41] lazyeval_0.2.0       RCurl_1.95-4.8       tibble_1.2           RSQLite_1.1         
[45] Formula_1.2-1        cluster_2.0.5        Matrix_1.2-8         data.table_1.9.8    
[49] assertthat_0.1       rpart_4.1-10         nnet_7.3-12 
ADD COMMENTlink written 2.8 years ago by Michael Love26k

Hi Michael I am still trying to figure out where the discrepancy occurs.  I had also been filtering using fpm and thought that was causing the error but I cannot reproduce the error when using the example DESeqDataSet.

ADD REPLYlink modified 2.7 years ago • written 2.7 years ago by Zach Roe10
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 16.09
Traffic: 398 users visited in the last hour