More significant genes with alpha 0.05 than 0.1 DESeq2
1
0
Entering edit mode
frene • 0
@frene-20306
Last seen 5.0 years ago

Hello everybody,

As I can read on the DESEQ2 manual, if we want an other p-adj diferent than 0.1 we have to select it by changing the alpha value. When I change it in the result formula:

alpha <- 0.05
res_0.05_filt <- results(dds1,contrast=c("condition","WT", "OE"),alpha=alpha)

Why I obtain more significant genes with p-adj 0.05 than 0.1? Why the mean count change in both cases? It is correct to do the results with 0.1 and them filter by 0.05?

here is the DESeq2 exit:

**> summary(res1)

out of 24732 with nonzero total read count
adjusted p-value < 0.1
LFC > 0 (up)       : 17, 0.069%
LFC < 0 (down)     : 24, 0.097%
outliers [1]       : 0, 0%
low counts [2]     : 6672, 27%
(mean count < 17)
[1] see 'cooksCutoff' argument of ?results
[2] see 'independentFiltering' argument of ?results

> table(res1$padj<=0.05)

FALSE  TRUE 
18050    10 
> summary(res_0.05_filt)

out of 24732 with nonzero total read count
adjusted p-value < 0.05
LFC > 0 (up)       : 8, 0.032%
LFC < 0 (down)     : 13, 0.053%
outliers [1]       : 0, 0%
low counts [2]     : 16679, 67%
(mean count < 239)
[1] see 'cooksCutoff' argument of ?results
[2] see 'independentFiltering' argument of ?results

> table(res_0.05_filt$padj<=0.05)

FALSE  TRUE 
 8032    21**

Here is the sessioninfo:

> sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows >= 8 x64 (build 9200)

Matrix products: default

locale:
[1] LC_COLLATE=Spanish_Spain.1252  LC_CTYPE=Spanish_Spain.1252    LC_MONETARY=Spanish_Spain.1252 LC_NUMERIC=C                  
[5] LC_TIME=Spanish_Spain.1252    

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

other attached packages:
 [1] ggbiplot_0.55               scales_1.0.0                plyr_1.8.4                  corrplot_0.84               RColorBrewer_1.1-2         
 [6] pheatmap_1.0.12             sva_3.30.1                  genefilter_1.64.0           mgcv_1.8-28                 nlme_3.1-137               
[11] gplots_3.0.1.1              ggplot2_3.1.0               limma_3.38.3                DESeq2_1.22.2               SummarizedExperiment_1.12.0
[16] DelayedArray_0.8.0          BiocParallel_1.16.6         matrixStats_0.54.0          Biobase_2.42.0              GenomicRanges_1.34.0       
[21] GenomeInfoDb_1.18.2         IRanges_2.16.0              S4Vectors_0.20.1            BiocGenerics_0.28.0        

loaded via a namespace (and not attached):
 [1] bit64_0.9-7            splines_3.5.1          gtools_3.8.1           Formula_1.2-3          latticeExtra_0.6-28    blob_1.1.1            
 [7] GenomeInfoDbData_1.2.0 pillar_1.3.1           RSQLite_2.1.1          backports_1.1.3        lattice_0.20-38        digest_0.6.18         
[13] XVector_0.22.0         checkmate_1.9.1        colorspace_1.4-1       htmltools_0.3.6        Matrix_1.2-17          XML_3.98-1.19         
[19] pkgconfig_2.0.2        zlibbioc_1.28.0        xtable_1.8-3           gdata_2.18.0           htmlTable_1.13.1       tibble_2.1.1          
[25] annotate_1.60.1        withr_2.1.2            nnet_7.3-12            lazyeval_0.2.2         survival_2.43-3        magrittr_1.5          
[31] crayon_1.3.4           memoise_1.1.0          MASS_7.3-51.1          foreign_0.8-71         tools_3.5.1            data.table_1.12.0     
[37] stringr_1.4.0          munsell_0.5.0          locfit_1.5-9.1         cluster_2.0.7-1        AnnotationDbi_1.44.0   compiler_3.5.1        
[43] caTools_1.17.1.2       rlang_0.3.2            RCurl_1.95-4.12        rstudioapi_0.10        htmlwidgets_1.3        labeling_0.3          
[49] bitops_1.0-6           base64enc_0.1-3        gtable_0.3.0           DBI_1.0.0              gridExtra_2.3          knitr_1.22            
[55] bit_1.1-14             Hmisc_4.2-0            KernSmooth_2.23-15     stringi_1.4.3          Rcpp_1.0.1             geneplotter_1.60.0    
[61] rpart_4.1-13           acepack_1.4.1          xfun_0.5              
>
deseq2 • 812 views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 21 minutes ago
United States

This is just what you would expect from the documentation. The alpha you provide to results() informs the independent filtering, and giving it the alpha you intend to use will give you the most rejections.

ADD COMMENT
0
Entering edit mode

Thanks Michel,

I saw a problem because with alpha 0.1 I obtained 10 genes and with 0.05, 21 genes, but it's true that the low counts rejection number of genes is greater with 0.05.

Another related question to de independent filtering. This is probably because of my limited knowledge , but If you turn it off with independentFiltering=FALSE, you get the p-adj values for all the genes but the values are higher than if you do with the independentFiltering=TRUE (default).
is there any way to make the independent Filtering with a basemean value defined by the user?

thanks a lot.

ADD REPLY
0
Entering edit mode

No it’s an automatic procedure (see paper or vignettte).

You can turn it off and filter the results table and then recompute adjusted pvalues with p.adjust().

ADD REPLY

Login before adding your answer.

Traffic: 656 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6