Entering edit mode
Dear friends,
I have been using diffbind to get RLE/TMM normalized counts. I have noticed that the RLE/TMM normalization is performed based on reads subtracted by control reads. How to turn off the subtraction? I have tried following by setting
bSubControl = F
, but it does not work.
library(DiffBind)
library(tidyverse)
h3k27ac_samples <- read.csv("h3k27ac.csv") #use sorted bam
h3k27ac <- dba(sampleSheet=h3k27ac_samples)
h3k27ac <- dba.blacklist(h3k27ac, blacklist =T , greylist = F)
h3k27ac <- dba.count(h3k27ac, bSubControl = F) #here the default score is library size norm
#get raw counts
h3k27ac <- dba.count(h3k27ac, peaks=NULL, score=DBA_SCORE_READS, bSubControl = F)
h3k27ac_count_raw <- dba.peakset(h3k27ac, bRetrieve=TRUE)
#get library size normalized reads
h3k27ac <- dba.normalize(h3k27ac) #default normalize = DBA_NORM_DEFAULT meaning normalize=DBA_NORM_LIB library=DBA_LIBSIZE_FULL background=FALSE
h3k27ac <- dba.count(h3k27ac, peaks = NULL, score = DBA_SCORE_NORMALIZED, bSubControl = F)
h3k27ac_count_lib_norm = dba.peakset(h3k27ac, bRetrieve=TRUE)
h3k27ac_lib_norm_info = dba.normalize(h3k27ac, bRetrieve = T)
# get RLE normalized counts
h3k27ac <- dba.normalize(h3k27ac, normalize = DBA_NORM_RLE) # normalize = DBA_NORM_RLE
h3k27ac <- dba.count(h3k27ac, peaks = NULL, score = DBA_SCORE_NORMALIZED, bSubControl = F)
h3k27ac_count_rle_norm = dba.peakset(h3k27ac, bRetrieve =TRUE)
h3k27ac_rle_norm_info <- dba.normalize(h3k27ac, normalize = DBA_NORM_RLE, bRetrieve = T)
# get TMM normalized counts
h3k27ac <- dba.normalize(h3k27ac, normalize = DBA_NORM_TMM) # normalize = DBA_NORM_TMM
h3k27ac <- dba.count(h3k27ac, peaks = NULL, score = DBA_SCORE_NORMALIZED, bSubControl = F)
h3k27ac_count_tmm_norm = dba.peakset(h3k27ac, bRetrieve=TRUE)
h3k27ac_tmm_norm_info <- dba.normalize(h3k27ac, normalize = DBA_NORM_TMM, bRetrieve = T)
None of h3k27ac_lib_norm_info, h3k27ac_rle_norm_info and h3k27ac_tmm_norm_info shows $control.subtract = F
for example
>h3k27ac_tmm_norm_info$control.subtract
[1] TRUE
Here is my sessionInfo
> sessionInfo()
R version 4.2.0 (2022-04-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 20.04.3 LTS
Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=C LC_ADDRESS=C LC_TELEPHONE=C LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] stats4 stats graphics grDevices utils datasets methods base
other attached packages:
[1] forcats_0.5.1 stringr_1.4.0 dplyr_1.0.9 purrr_0.3.4 readr_2.1.2 tidyr_1.2.0
[7] tibble_3.1.7 ggplot2_3.3.6 tidyverse_1.3.1 DiffBind_3.6.1 SummarizedExperiment_1.26.1 Biobase_2.56.0
[13] MatrixGenerics_1.8.0 matrixStats_0.62.0 GenomicRanges_1.48.0 GenomeInfoDb_1.32.2 IRanges_2.30.0 S4Vectors_0.34.0
[19] BiocGenerics_0.42.0
loaded via a namespace (and not attached):
[1] amap_0.8-18 colorspace_2.0-3 rjson_0.2.21 hwriter_1.3.2.1 ellipsis_0.3.2 XVector_0.36.0
[7] fs_1.5.2 rstudioapi_0.13 bit64_4.0.5 ggrepel_0.9.1 AnnotationDbi_1.58.0 fansi_1.0.3
[13] mvtnorm_1.1-3 apeglm_1.18.0 lubridate_1.8.0 xml2_1.3.3 splines_4.2.0 cachem_1.0.6
[19] geneplotter_1.74.0 jsonlite_1.8.0 Rsamtools_2.12.0 annotate_1.74.0 broom_0.8.0 ashr_2.2-54
[25] dbplyr_2.1.1 png_0.1-7 GreyListChIP_1.28.1 compiler_4.2.0 httr_1.4.3 backports_1.4.1
[31] assertthat_0.2.1 Matrix_1.4-1 fastmap_1.1.0 limma_3.52.1 cli_3.3.0 htmltools_0.5.2
[37] tools_4.2.0 coda_0.19-4 gtable_0.3.0 glue_1.6.2 GenomeInfoDbData_1.2.8 systemPipeR_2.2.2
[43] ShortRead_1.54.0 Rcpp_1.0.8.3 bbmle_1.0.25 cellranger_1.1.0 vctrs_0.4.1 Biostrings_2.64.0
[49] rtracklayer_1.56.0 rvest_1.0.2 lifecycle_1.0.1 irlba_2.3.5 restfulr_0.0.13 gtools_3.9.2.1
[55] XML_3.99-0.9 edgeR_3.38.1 zlibbioc_1.42.0 MASS_7.3-57 scales_1.2.0 BSgenome_1.64.0
[61] hms_1.1.1 parallel_4.2.0 RColorBrewer_1.1-3 yaml_2.3.5 memoise_2.0.1 emdbook_1.3.12
[67] bdsmatrix_1.3-4 RSQLite_2.2.14 latticeExtra_0.6-29 stringi_1.7.6 SQUAREM_2021.1 genefilter_1.78.0
[73] BiocIO_1.6.0 caTools_1.18.2 BiocParallel_1.30.2 truncnorm_1.0-8 rlang_1.0.2 pkgconfig_2.0.3
[79] bitops_1.0-7 lattice_0.20-45 invgamma_1.1 GenomicAlignments_1.32.0 htmlwidgets_1.5.4 bit_4.0.4
[85] tidyselect_1.1.2 DESeq2_1.36.0 plyr_1.8.7 magrittr_2.0.3 R6_2.5.1 gplots_3.1.3
[91] generics_0.1.2 DelayedArray_0.22.0 DBI_1.1.2 pillar_1.7.0 haven_2.5.0 withr_2.5.0
[97] survival_3.3-1 KEGGREST_1.36.0 RCurl_1.98-1.6 mixsqp_0.3-43 modelr_0.1.8 crayon_1.5.1
[103] KernSmooth_2.23-20 utf8_1.2.2 tzdb_0.3.0 jpeg_0.1-9 locfit_1.5-9.5 grid_4.2.0
[109] readxl_1.4.0 blob_1.2.3 reprex_2.0.1 digest_0.6.29 xtable_1.8-4 numDeriv_2016.8-1.1
[115] munsell_0.5.0