Question

Error when using lfcShrink and contrast argument in DESeq2

0

Entering edit mode

Oskar ▴ 50

@oskar-13385

Last seen 5.8 years ago

Hi everyone -

I'm running a Differential Expression analysis via DESeq2 (DESeq2_1.16.1) .

I have three levels in the conditions (SNP_A, SNP_B, and control), and I set as reference the "control" level. Each condition has three biological replicates. Then I ran this:

dds <- DESeqDataSet(se, design = ~ condition)
dds <- DESeq(dds)

From the dds obtained above, I ran pairwise comparison between conditions as following:

ACres <- results(dds, contrast=c("condition", "SNP_A", "control"))
BCres <- results(dds, contrast=c("condition", "SNP_B", "control"))
BAres <- results(dds, contrast=c("condition", "SNP_B", "SNP_A"))

So far, all work well. However, when I try to run the lfcShrink function for the same contrast, only work for the first two (the ones comparing each SNP with the control), but not for the contrasts between SNPs. Here the code I used:

ACres.lfc <- lfcShrink(dds, contrast=c("condition", "SNP_A", "control"), res = ACres)
BCres.lfc <- lfcShrink(dds, contrast=c("condition", "SNP_B", "control"), res = BCres)
BAres.lfc <- lfcShrink(dds, contrast=c("condition", "SNP_B", "SNP_A"), res = BAres)

I got the following error:

Error in cleanContrast(object, contrast, expanded = isExpanded, listValues = listValues,  : 
  conditionSNP_B and conditionSNP_A are expected to be in resultsNames(object)

When I check the resultsNames(dds) I got only

[1] "Intercept" "condition_SNP_A_vs_control" "condition_SNP_B_vs_control"

I don't understand why I can not get all possible contrasts when using the lfcShrink function as when using only the results function. How can I run this other contrast? Should I set another level as a reference and redoing all again?

Many thanks in advance in this matter.

deseq2 lfcshrink • 2.3k views

ADD COMMENT • link 7.3 years ago Oskar ▴ 50

0

Entering edit mode

Hi Michael -

I used "dots" instead underscores and it works. Thank you very much indeed.

ADD REPLY • link 7.3 years ago Oskar ▴ 50

0

Entering edit mode

hi Oskar,

If you have a chance, can you try my example below? I'm confused as to why I don't see the same error. (And would like to fix the bug if it exists for underscores in factor levels and lfcShrink usage). Was everything run with DESeq2 v1.16.1 from Bioconductor?

ADD REPLY • link 7.3 years ago Michael Love 42k

score 1 · Answer 1 · 2017-08-18

1

Entering edit mode

Michael Love 42k

@mikelove

Last seen 20 hours ago

United States

Can you try using factor levels without underscore "_". I can investigate as well but that might fix it.

ADD COMMENT • link 7.3 years ago Michael Love 42k

0

Entering edit mode

I can't reproduce this bug on my end using underscores. Does the following example work for you?

dds <- makeExampleDESeqDataSet()
dds$condition <- factor(rep(c("control","SNP_A","SNP_B"),each=4))
dds <- DESeq(dds)
res <- results(dds, contrast=c("condition","SNP_B","SNP_A"))
res <- lfcShrink(dds, contrast=c("condition","SNP_B","SNP_A"), res=res)

> sessionInfo()
R version 3.4.1 (2017-06-30)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 17.04

Matrix products: default
BLAS: /usr/lib/libblas/libblas.so.3.7.0
LAPACK: /usr/lib/lapack/liblapack.so.3.7.0

locale:
[1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C             
[3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8   
[5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8  
[7] LC_PAPER=en_US.UTF-8       LC_NAME=C                
[9] LC_ADDRESS=C               LC_TELEPHONE=C           
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

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

other attached packages:
[1] DESeq2_1.16.1              SummarizedExperiment_1.6.3
[3] DelayedArray_0.2.7         matrixStats_0.52.2       
[5] Biobase_2.36.2             GenomicRanges_1.28.3     
[7] GenomeInfoDb_1.12.2        IRanges_2.10.2           
[9] S4Vectors_0.14.3           BiocGenerics_0.22.0      
[11] rmarkdown_1.6              magrittr_1.5             
[13] testthat_1.0.2             devtools_1.13.2           

loaded via a namespace (and not attached):
[1] locfit_1.5-9.1          Rcpp_0.12.11            lattice_0.20-35       
[4] rprojroot_1.2           digest_0.6.12           R6_2.2.2              
[7] plyr_1.8.4              backports_1.1.0         acepack_1.4.1         
[10] RSQLite_2.0             evaluate_0.10.1         ggplot2_2.2.1         
[13] zlibbioc_1.22.0         rlang_0.1.1             lazyeval_0.2.0        
[16] data.table_1.10.4       annotate_1.54.0         blob_1.1.0            
[19] rpart_4.1-11            Matrix_1.2-10           checkmate_1.8.2       
[22] splines_3.4.1           BiocParallel_1.10.1     geneplotter_1.54.0    
[25] stringr_1.2.0           foreign_0.8-69          htmlwidgets_0.8       
[28] bit_1.1-12              RCurl_1.95-4.8          munsell_0.4.3         
[31] compiler_3.4.1          base64enc_0.1-3         htmltools_0.3.6       
[34] nnet_7.3-12             tibble_1.3.3            gridExtra_2.2.1       
[37] htmlTable_1.9           GenomeInfoDbData_0.99.0 Hmisc_4.0-3           
[40] XML_3.98-1.9            crayon_1.3.2            withr_1.0.2           
[43] bitops_1.0-6            grid_3.4.1              xtable_1.8-2          
[46] gtable_0.2.0            DBI_0.7                 scales_0.4.1          
[49] stringi_1.1.5           XVector_0.16.0          genefilter_1.58.1     
[52] latticeExtra_0.6-28     Formula_1.2-1           RColorBrewer_1.1-2    
[55] tools_3.4.1             bit64_0.9-7             survival_2.41-3       
[58] AnnotationDbi_1.38.1    colorspace_1.3-2        cluster_2.0.6         
[61] memoise_1.1.0           knitr_1.16

ADD REPLY • link 7.3 years ago Michael Love 42k