Question: DiffBind dba.count error pv$merge(bed): attempt to apply non-function 0 6 months ago by chadnewton0 wrote: I am attempting to perform differential peak analysis from ATAC-seq data using DiffBind package and am running into problems with the dba.count function. I have 2 conditions (rapid vs slow) with 3 samples each for a total of 6 samples. I want to create a consensus peakset for both rapid and slow condition. I would like to identify punctate peaks, so I set the summits = 200. cn <- dba(sampleSheet = "cn_atac.csv", peakFormat = "bed", minOverlap = 1, config = data.frame(RunParallel=TRUE, reportInit='DBA', DataType=DBA_DATA_GRANGES, minQCth=15, fragmentSize=125, bCorPlot=FALSE, th=0.05, bUsePval=FALSE), bRemoveM = TRUE, bRemoveRandom = TRUE, filter = 20) cn.consensus = dba.peakset(cn, consensus = DBA_CONDITION, minOverlap = 1) peaks = dba.peakset(cn.consensus, cn.consensus$masks$Consensus, bRetrieve = TRUE) cn.cons.count = dba.count(cn, peaks = peaks, summits = 200, bRemoveDuplicates = FALSE, bScaleControl = TRUE, bParallel = FALSE)  generates error:  Sample: atac_files/c1.bam125 [W::bam_hdr_read] bgzf_check_EOF: No error [W::bam_hdr_read] bgzf_check_EOF: No error Sample: atac_files/c2.bam125 Sample: atac_files/c3.bam125 [W::bam_hdr_read] bgzf_check_EOF: No error [W::bam_hdr_read] bgzf_check_EOF: No error Sample: atac_files/c4.bam125 Sample: atac_files/c5.bam125 [W::bam_hdr_read] bgzf_check_EOF: No error [W::bam_hdr_read] bgzf_check_EOF: No error Sample: atac_files/c6.bam125 Re-centering peaks... Error in pv$merge(bed) : attempt to apply non-function
In addition: Warning message:
In dba.multicore.init(DBA$config) : Parallel execution unavailable: executing serially. > traceback() 3: pv.Recenter(pv, summits, (numpeaks - numAdded + 1):numpeaks, called) 2: pv.counts(DBA, peaks = peaks, minOverlap = minOverlap, defaultScore = score, bLog = bLog, insertLength = fragmentSize, bOnlyCounts = T, bCalledMasks = TRUE, minMaxval = filter, bParallel = bParallel, bUseLast = bUseLast, bWithoutDupes = bRemoveDuplicates, bScaleControl = bScaleControl, filterFun = filterFun, bLowMem = bUseSummarizeOverlaps, readFormat = readFormat, summits = summits, minMappingQuality = mapQCth) 1: dba.count(cn, peaks = peaks, summits = 200, bRemoveDuplicates = FALSE, bScaleControl = TRUE, bParallel = FALSE)  . > sessionInfo() R version 3.6.0 (2019-04-26) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 17134) Matrix products: default locale: [1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 [3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C [5] LC_TIME=English_United States.1252 attached base packages: [1] parallel stats4 stats graphics grDevices utils datasets methods base other attached packages: [1] DiffBind_2.12.0 SummarizedExperiment_1.14.0 DelayedArray_0.10.0 [4] BiocParallel_1.17.18 matrixStats_0.54.0 Biobase_2.44.0 [7] GenomicRanges_1.36.0 GenomeInfoDb_1.20.0 IRanges_2.18.0 [10] S4Vectors_0.22.0 BiocGenerics_0.30.0 loaded via a namespace (and not attached): [1] amap_0.8-17 colorspace_1.4-1 rjson_0.2.20 hwriter_1.3.2 [5] htmlTable_1.13.1 XVector_0.24.0 base64enc_0.1-3 rstudioapi_0.10 [9] ggrepel_0.8.1 bit64_0.9-7 AnnotationDbi_1.46.0 splines_3.6.0 [13] geneplotter_1.62.0 knitr_1.23 Formula_1.2-3 Rsamtools_2.0.0 [17] annotate_1.62.0 cluster_2.0.8 GO.db_3.8.2 pheatmap_1.0.12 [21] graph_1.62.0 BiocManager_1.30.4 compiler_3.6.0 httr_1.4.0 [25] GOstats_2.50.0 backports_1.1.4 assertthat_0.2.1 Matrix_1.2-17 [29] lazyeval_0.2.2 limma_3.40.2 htmltools_0.3.6 acepack_1.4.1 [33] prettyunits_1.0.2 tools_3.6.0 gtable_0.3.0 glue_1.3.1 [37] GenomeInfoDbData_1.2.1 Category_2.50.0 systemPipeR_1.18.0 dplyr_0.8.1 [41] batchtools_0.9.11 rappdirs_0.3.1 ShortRead_1.42.0 Rcpp_1.0.1 [45] Biostrings_2.52.0 gdata_2.18.0 rtracklayer_1.44.0 xfun_0.7 [49] stringr_1.4.0 gtools_3.8.1 XML_3.98-1.19 edgeR_3.26.3 [53] zlibbioc_1.30.0 scales_1.0.0 BSgenome_1.52.0 VariantAnnotation_1.30.1 [57] hms_0.4.2 RBGL_1.60.0 RColorBrewer_1.1-2 yaml_2.2.0 [61] gridExtra_2.3 memoise_1.1.0 ggplot2_3.1.1 biomaRt_2.40.0 [65] rpart_4.1-15 latticeExtra_0.6-28 stringi_1.4.3 RSQLite_2.1.1 [69] genefilter_1.66.0 checkmate_1.9.3 GenomicFeatures_1.36.0 caTools_1.17.1.2 [73] rlang_0.3.4 pkgconfig_2.0.2 bitops_1.0-6 lattice_0.20-38 [77] purrr_0.3.2 labeling_0.3 htmlwidgets_1.3 GenomicAlignments_1.20.0 [81] bit_1.1-14 tidyselect_0.2.5 GSEABase_1.46.0 AnnotationForge_1.26.0 [85] plyr_1.8.4 magrittr_1.5 DESeq2_1.24.0 R6_2.4.0 [89] gplots_3.0.1.1 Hmisc_4.2-0 base64url_1.4 DBI_1.0.0 [93] pillar_1.4.0 foreign_0.8-71 withr_2.1.2 survival_2.44-1.1 [97] RCurl_1.95-4.12 nnet_7.3-12 tibble_2.1.1 crayon_1.3.4 [101] KernSmooth_2.23-15 progress_1.2.2 locfit_1.5-9.1 grid_3.6.0 [105] data.table_1.12.2 blob_1.1.1 Rgraphviz_2.28.0 digest_0.6.19 [109] xtable_1.8-4 brew_1.0-6 munsell_0.5.0  Much appreciated! Chad diffbind • 385 views ADD COMMENTlink modified 3 months ago by Rory Stark3.0k • written 6 months ago by chadnewton0 Did you figure this out? I'm having much the same issue. There's no discernible difference in our code so I wont dump that on here, but I have included my session info below if anyone else needs it. Another similarity is that I am trying to use this for ATACseq. With that I do not having any inputs. Could this contribute? > sessionInfo() R version 3.6.0 (2019-04-26) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Ubuntu 18.04.2 LTS Matrix products: default BLAS: /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3 LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so Random number generation: RNG: Mersenne-Twister Normal: Inversion Sample: Rounding 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 [6] LC_MESSAGES=en_US.UTF-8 LC_PAPER=en_US.UTF-8 LC_NAME=C 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 utils datasets methods base other attached packages: [1] DiffBind_2.12.0 SummarizedExperiment_1.14.0 DelayedArray_0.10.0 BiocParallel_1.18.0 [5] matrixStats_0.54.0 Biobase_2.44.0 GenomicRanges_1.36.0 GenomeInfoDb_1.20.0 [9] IRanges_2.18.1 S4Vectors_0.22.0 BiocGenerics_0.30.0 BiFET_1.4.0 loaded via a namespace (and not attached): [1] Category_2.50.0 bitops_1.0-6 bit64_0.9-7 RColorBrewer_1.1-2 progress_1.2.2 [6] httr_1.4.0 Rgraphviz_2.28.0 backports_1.1.4 tools_3.6.0 R6_2.4.0 [11] KernSmooth_2.23-15 DBI_1.0.0 lazyeval_0.2.2 colorspace_1.4-1 withr_2.1.2 [16] tidyselect_0.2.5 prettyunits_1.0.2 bit_1.1-14 compiler_3.6.0 graph_1.62.0 [21] rtracklayer_1.44.0 checkmate_1.9.3 caTools_1.17.1.2 scales_1.0.0 genefilter_1.66.0 [26] RBGL_1.60.0 rappdirs_0.3.1 stringr_1.4.0 digest_0.6.19 Rsamtools_2.0.0 [31] AnnotationForge_1.26.0 XVector_0.24.0 pkgconfig_2.0.2 BSgenome_1.52.0 limma_3.40.2 [36] rlang_0.3.4 rstudioapi_0.10 RSQLite_2.1.1 GOstats_2.50.0 hwriter_1.3.2 [41] gtools_3.8.1 dplyr_0.8.1 VariantAnnotation_1.30.1 RCurl_1.95-4.12 magrittr_1.5 [46] GO.db_3.8.2 GenomeInfoDbData_1.2.1 Matrix_1.2-17 Rcpp_1.0.1 munsell_0.5.0 [51] stringi_1.4.3 yaml_2.2.0 edgeR_3.26.4 zlibbioc_1.30.0 gplots_3.0.1.1 [56] plyr_1.8.4 grid_3.6.0 blob_1.1.1 ggrepel_0.8.1 gdata_2.18.0 [61] crayon_1.3.4 lattice_0.20-38 splines_3.6.0 Biostrings_2.52.0 GenomicFeatures_1.36.1 [66] annotate_1.62.0 poibin_1.3 hms_0.4.2 batchtools_0.9.11 locfit_1.5-9.1 [71] knitr_1.23 pillar_1.4.1 rjson_0.2.20 systemPipeR_1.18.1 base64url_1.4 [76] biomaRt_2.40.0 XML_3.98-1.20 glue_1.3.1 ShortRead_1.42.0 latticeExtra_0.6-28 [81] data.table_1.12.2 gtable_0.3.0 purrr_0.3.2 amap_0.8-17 assertthat_0.2.1 [86] ggplot2_3.1.1 xfun_0.7 xtable_1.8-4 survival_2.44-1.1 pheatmap_1.0.12 [91] tibble_2.1.3 GenomicAlignments_1.20.0 AnnotationDbi_1.46.0 memoise_1.1.0 brew_1.0-6 [96] GSEABase_1.46.0  ADD REPLYlink written 6 months ago by matthew.paul.20060 Answer: DiffBind dba.count error pv$merge(bed): attempt to apply non-function
0
3 months ago by
Rory Stark3.0k
CRUK, Cambridge, UK
Rory Stark3.0k wrote:

I'm suspicous that some of the files generate a [W::bam_hdr_read] bgzf_check_EOF: No error warning and some do not.

Can you try calling dba.count() a) without setting summits or b) setting summits=TRUE? Does it work in either of these cases? If it works with summits=TRUE, you can run dba.count() again with peaks=NULL and summits=200, and see if that works.

If it works with summits=TRUE, but fails when you call it again with peaks=NULL and summits=200, please send me a link to the DBA object after calling with with summits=TRUE so I can have a look.

-Rory