dba.count() error when re-centering peaks
1
0
Entering edit mode
Mark • 0
@mark-24951
Last seen 8 months ago
United States

Hi there,

I got the error message below when using dba.count() with consensus peaks. Can you help to solve the problem? Thanks!

Re-centering peaks...
Error in heights * sapply(called, function(x) x) : non-conformable arrays

Here is more information. My code:


data(tamoxifen_peaks)
tamoxifen_consensus <-
  dba.peakset(tamoxifen,
              consensus = c(DBA_TISSUE, DBA_CONDITION),
              minOverlap = 0.66)

tamoxifen_consensus <- dba(tamoxifen_consensus,
                           mask = tamoxifen_consensus$masks$Consensus,
                           minOverlap = 1)

consensus_peaks <- dba.peakset(tamoxifen_consensus, bRetrieve = TRUE)

setwd(vignette_data_dir)
tamoxifen <- dba.count(tamoxifen, peaks=consensus_peaks)

Here is the output and error message:

Computing summits...
Sample: reads/Chr18_BT474_ER_1.bam125 
Sample: reads/Chr18_BT474_ER_2.bam125 
Sample: reads/Chr18_MCF7_ER_1.bam125 
Sample: reads/Chr18_MCF7_ER_2.bam125 
Sample: reads/Chr18_MCF7_ER_3.bam125 
Sample: reads/Chr18_T47D_ER_1.bam125 
Sample: reads/Chr18_T47D_ER_2.bam125 
Sample: reads/Chr18_TAMR_ER_1.bam125 
Sample: reads/Chr18_TAMR_ER_2.bam125 
Sample: reads/Chr18_ZR75_ER_1.bam125 
Sample: reads/Chr18_ZR75_ER_2.bam125 
Sample: reads/Chr18_BT474_input.bam125 
Sample: reads/Chr18_MCF7_input.bam125 
Sample: reads/Chr18_T47D_input.bam125 
Sample: reads/Chr18_TAMR_input.bam125 
Sample: reads/Chr18_ZR75_input.bam125 
Re-centering peaks...
Error in heights * sapply(called, function(x) x) : non-conformable arrays

traceback() returns:

3: pv.Recenter(pv, summits, (numpeaks - numAdded + 1):numpeaks, 
       pv$called)
2: pv.counts(DBA, peaks = peaks, minOverlap = minOverlap, defaultScore = score, 
       bLog = bLog, insertLength = fragmentSize, bOnlyCounts = TRUE, 
       bCalledMasks = TRUE, minMaxval = filter, bParallel = bParallel, 
       bUseLast = bUseLast, bWithoutDupes = bRemoveDuplicates, bScaleControl = bScaleControl, 
       filterFun = filterFun, bLowMem = bUseSummarizeOverlaps, readFormat = readFormat, 
       summits = summits, minMappingQuality = mapQCth, minCount = minCount, 
       bSubControl = bSubControl, maxGap = maxGap)
1: dba.count(tamoxifen, peaks = consensus_peaks)

Session information:

R version 4.1.2 (2021-11-01)
Platform: x86_64-apple-darwin13.4.0 (64-bit)
Running under: macOS Big Sur 10.16

Matrix products: default
LAPACK: /usr/local/anaconda3/envs/bioc3.14-py39/lib/libopenblasp-r0.3.18.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

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

other attached packages:
 [1] forcats_0.5.1               stringr_1.4.0              
 [3] dplyr_1.0.7                 purrr_0.3.4                
 [5] readr_2.1.1                 tidyr_1.1.4                
 [7] tibble_3.1.6                ggplot2_3.3.5              
 [9] tidyverse_1.3.1             DiffBind_3.4.7             
[11] csaw_1.28.0                 SummarizedExperiment_1.24.0
[13] Biobase_2.54.0              MatrixGenerics_1.6.0       
[15] matrixStats_0.61.0          GenomicRanges_1.46.1       
[17] GenomeInfoDb_1.30.0         IRanges_2.28.0             
[19] S4Vectors_0.32.3            BiocGenerics_0.40.0        
[21] shiny_1.7.1                

loaded via a namespace (and not attached):
  [1] readxl_1.3.1             backports_1.4.1         
  [3] servr_0.24               plyr_1.8.6              
  [5] BiocParallel_1.28.3      usethis_2.1.5           
  [7] amap_0.8-18              digest_0.6.29           
  [9] invgamma_1.1             htmltools_0.5.2         
 [11] SQUAREM_2021.1           fansi_1.0.2             
 [13] magrittr_2.0.2           memoise_2.0.1           
 [15] BSgenome_1.62.0          tzdb_0.2.0              
 [17] limma_3.50.0             remotes_2.4.2           
 [19] Biostrings_2.62.0        modelr_0.1.8            
 [21] systemPipeR_2.0.5        bdsmatrix_1.3-4         
 [23] prettyunits_1.1.1        jpeg_0.1-9              
 [25] colorspace_2.0-2         rvest_1.0.2             
 [27] apeglm_1.16.0            ggrepel_0.9.1           
 [29] haven_2.4.3              xfun_0.29               
 [31] callr_3.7.0              crayon_1.4.2            
 [33] RCurl_1.98-1.5           jsonlite_1.7.3          
 [35] glue_1.6.1               gtable_0.3.0            
 [37] zlibbioc_1.40.0          XVector_0.34.0          
 [39] DelayedArray_0.20.0      pkgbuild_1.3.1          
 [41] scales_1.1.1             mvtnorm_1.1-3           
 [43] DBI_1.1.2                edgeR_3.36.0            
 [45] miniUI_0.1.1.1           Rcpp_1.0.8              
 [47] xtable_1.8-4             emdbook_1.3.12          
 [49] truncnorm_1.0-8          httr_1.4.2              
 [51] metapod_1.2.0            htmlwidgets_1.5.4       
 [53] gplots_3.1.1             RColorBrewer_1.1-2      
 [55] ellipsis_0.3.2           pkgconfig_2.0.3         
 [57] XML_3.99-0.8             dbplyr_2.1.1            
 [59] sass_0.4.0               locfit_1.5-9.4          
 [61] utf8_1.2.2               tidyselect_1.1.1        
 [63] rlang_1.0.0              later_1.3.0             
 [65] cellranger_1.1.0         munsell_0.5.0           
 [67] tools_4.1.2              cachem_1.0.6            
 [69] cli_3.1.1                generics_0.1.1          
 [71] broom_0.7.12             devtools_2.4.3          
 [73] evaluate_0.14            fastmap_1.1.0           
 [75] yaml_2.2.2               processx_3.5.2          
 [77] knitr_1.37               fs_1.5.2                
 [79] caTools_1.18.2           mime_0.12               
 [81] xml2_1.3.3               brio_1.1.3              
 [83] compiler_4.1.2           rstudioapi_0.13         
 [85] png_0.1-7                testthat_3.1.2          
 [87] reprex_2.0.1             bslib_0.3.1             
 [89] stringi_1.7.6            ps_1.6.0                
 [91] blogdown_1.7             desc_1.4.0              
 [93] lattice_0.20-45          Matrix_1.4-0            
 [95] vctrs_0.3.8              pillar_1.6.5            
 [97] lifecycle_1.0.1          BiocManager_1.30.16     
 [99] jquerylib_0.1.4          bitops_1.0-7            
[101] irlba_2.3.5              httpuv_1.6.5            
[103] rtracklayer_1.54.0       BiocIO_1.4.0            
[105] R6_2.5.1                 latticeExtra_0.6-29     
[107] hwriter_1.3.2            bookdown_0.24           
[109] promises_1.2.0.1         ShortRead_1.52.0        
[111] KernSmooth_2.23-20       sessioninfo_1.2.2       
[113] MASS_7.3-55              gtools_3.9.2            
[115] assertthat_0.2.1         pkgload_1.2.4           
[117] rjson_0.2.21             rprojroot_2.0.2         
[119] withr_2.4.3              GenomicAlignments_1.30.0
[121] Rsamtools_2.10.0         GenomeInfoDbData_1.2.7  
[123] hms_1.1.1                parallel_4.1.2          
[125] grid_4.1.2               coda_0.19-4             
[127] rmarkdown_2.11           GreyListChIP_1.26.0     
[129] ashr_2.2-47              mixsqp_0.3-43           
[131] bbmle_1.0.24             lubridate_1.8.0         
[133] numDeriv_2016.8-1.1      restfulr_0.0.13
DiffBind • 718 views
ADD COMMENT
0
Entering edit mode
Rory Stark ★ 4.5k
@rory-stark-5741
Last seen 41 minutes ago
CRUK, Cambridge, UK

Thank you for a clear and complete report of your issue!

I have checked in a fix that will be available in the next few days as DiffBind_3.4.8.

ADD COMMENT
0
Entering edit mode

Thank you, Rory!

ADD REPLY
0
Entering edit mode

Hi Rory,

I am currently encountering the same error. I have copied my code below:

Tcf4_CRISPRa_ATAC <- dba(sampleSheet = "/Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/NJR0019_DiffBind.csv", minOverlap = 0)
Tcf4_CRISPRa_ATAC <- dba.count(Tcf4_CRISPRa_ATAC, summits = 250, bUseSummarizeOverlaps = FALSE, score = DBA_SCORE_TMM_READS_FULL_CPM)

The CSV file whose path is shown in the first line contains path info for all my BAM and BED files (I can provide the info from that file if it would be helpful). After running for a couple hours with no message, I received the following:

Re-centering peaks...
Error in heights * sapply(called, function(x) x) : non-conformable arrays

I have run this code as recently as 3 or 4 months ago with no issues. Using traceback() produced exactly the same message Mark saw/provided above.

Below is my session info:

sessionInfo()
R version 4.1.0 (2021-05-18)
Platform: x86_64-apple-darwin17.0 (64-bit)
Running under: macOS Big Sur 11.3.1

Matrix products: default
LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

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

other attached packages:
 [1] viridis_0.6.2               viridisLite_0.4.0           pheatmap_1.0.12             ggrepel_0.9.1              
 [5] ggplot2_3.3.6               sp_1.5-0                    SeuratObject_4.1.0          Seurat_4.1.1               
 [9] reshape2_1.4.4              rtracklayer_1.54.0          DiffBind_3.4.11             SummarizedExperiment_1.24.0
[13] Biobase_2.54.0              MatrixGenerics_1.6.0        matrixStats_0.61.0          GenomicRanges_1.46.1       
[17] GenomeInfoDb_1.30.1         IRanges_2.28.0              S4Vectors_0.32.4            BiocGenerics_0.40.0        

loaded via a namespace (and not attached):
  [1] plyr_1.8.7               igraph_1.2.11            lazyeval_0.2.2           splines_4.1.0            BiocParallel_1.28.3     
  [6] listenv_0.8.0            scattermore_0.8          amap_0.8-18              digest_0.6.29            invgamma_1.1            
 [11] htmltools_0.5.2          SQUAREM_2021.1           fansi_1.0.3              magrittr_2.0.2           BSgenome_1.62.0         
 [16] tensor_1.5               cluster_2.1.2            ROCR_1.0-11              limma_3.50.1             globals_0.14.0          
 [21] Biostrings_2.62.0        systemPipeR_2.0.5        spatstat.sparse_2.1-0    bdsmatrix_1.3-4          jpeg_0.1-9              
 [26] colorspace_2.0-3         apeglm_1.16.0            dplyr_1.0.8              crayon_1.5.0             RCurl_1.98-1.6          
 [31] jsonlite_1.8.0           spatstat.data_2.1-2      progressr_0.10.0         survival_3.3-1           zoo_1.8-9               
 [36] glue_1.6.2               polyclip_1.10-0          gtable_0.3.0             zlibbioc_1.40.0          XVector_0.34.0          
 [41] leiden_0.3.9             DelayedArray_0.20.0      future.apply_1.8.1       abind_1.4-5              scales_1.1.1            
 [46] mvtnorm_1.1-3            DBI_1.1.2                spatstat.random_2.1-0    miniUI_0.1.1.1           Rcpp_1.0.8.3            
 [51] xtable_1.8-4             emdbook_1.3.12           spatstat.core_2.4-0      reticulate_1.24          truncnorm_1.0-8         
 [56] htmlwidgets_1.5.4        httr_1.4.2               gplots_3.1.1             RColorBrewer_1.1-2       ellipsis_0.3.2          
 [61] ica_1.0-2                pkgconfig_2.0.3          XML_3.99-0.9             uwot_0.1.11              deldir_1.0-6            
 [66] locfit_1.5-9.5           utf8_1.2.2               tidyselect_1.1.2         rlang_1.0.2              later_1.3.0             
 [71] munsell_0.5.0            tools_4.1.0              cli_3.2.0                generics_0.1.2           ggridges_0.5.3          
 [76] stringr_1.4.0            fastmap_1.1.0            goftest_1.2-3            yaml_2.3.5               fitdistrplus_1.1-8      
 [81] caTools_1.18.2           purrr_0.3.4              RANN_2.6.1               nlme_3.1-155             pbapply_1.5-0           
 [86] future_1.24.0            mime_0.12                compiler_4.1.0           rstudioapi_0.13          plotly_4.10.0           
 [91] png_0.1-7                spatstat.utils_2.3-0     tibble_3.1.6             stringi_1.7.6            rgeos_0.5-9             
 [96] lattice_0.20-45          Matrix_1.4-1             vctrs_0.3.8              pillar_1.7.0             lifecycle_1.0.1         
[101] spatstat.geom_2.3-2      lmtest_0.9-40            RcppAnnoy_0.0.19         data.table_1.14.2        cowplot_1.1.1           
[106] bitops_1.0-7             irlba_2.3.5              httpuv_1.6.5             patchwork_1.1.1          R6_2.5.1                
[111] BiocIO_1.4.0             latticeExtra_0.6-29      hwriter_1.3.2            promises_1.2.0.1         ShortRead_1.52.0        
[116] gridExtra_2.3            KernSmooth_2.23-20       parallelly_1.30.0        codetools_0.2-18         MASS_7.3-56             
[121] gtools_3.9.2             assertthat_0.2.1         rjson_0.2.21             withr_2.5.0              GenomicAlignments_1.30.0
[126] sctransform_0.3.3        Rsamtools_2.10.0         GenomeInfoDbData_1.2.7   mgcv_1.8-39              parallel_4.1.0          
[131] rpart_4.1.16             grid_4.1.0               tidyr_1.2.0              coda_0.19-4              GreyListChIP_1.26.0     
[136] ashr_2.2-54              Rtsne_0.15               mixsqp_0.3-43            bbmle_1.0.24             numDeriv_2016.8-1.1     
[141] shiny_1.7.1              restfulr_0.0.13

I appreciate any help you can provide!

Nathaniel

ADD REPLY
0
Entering edit mode

Given that it runs for hours, I assume you have some combination of a large number of samples, peaks, and/or reads, is that correct?

The two things I recommend here are a) running this with the currently released version and b) running dba.count() with bParallel = FALSE. This may cause it to take an even longer time to run (although there are performance improvements in the current version over the one you are using).

One possibility is that DiffBind is taking too many cores for parallel processing (not enough memory to support all the parallel operations). So you could try setting Tcf4_CRISPRa_ATAC$config$cores <- 4 (or another number that works) to get some benefits of parallelization without overloading memory.

Another possibility: I recently had a report of someone having this problem. Then he told me how he fixed it:

I made a mistake while setting up the metadata table and had a duplicated BAM file in the bamReads column. Once I have replaced the duplicated file with the correct file for that particular sample, the analysis ran without any errors.

I am in the process of addressing this in an upcoming fix, but in the meantime you may want to check that your bamReads are all unique. If you need to use the same bamRead file for more than one sample, I have some workarounds for that.

If these do not solve your problem, let me know and we can figure out how to take a deeper look.

ADD REPLY
0
Entering edit mode

Hi Rory,

My apologies for the delay. I tried all of the measures you suggested, and I am still receiving the same error. Specifically:

  1. I updated to (what I believe is) the most recent version of Diffbind (v3.6.2; see Session Info below for details of all program versions)
  2. I and one of my colleagues double-checked the file paths in my metadata table, and all were correct
  3. I tried running with 4 CPUs rather than 8
  4. I added bParallel = FALSE to the dba.count function

Specifically for point (4), I ran:

> Tcf4_CRISPRa_ATAC <- dba.count(Tcf4_CRISPRa_ATAC, summits = 250, bUseSummarizeOverlaps = FALSE, score = DBA_SCORE_TMM_READS_FULL_CPM, bParallel = FALSE)

with the Tcf4_CRISPRa_ATAC object created as it was in my initial post. It may be worth noting that this did run faster than my initial attempt, even with the lack of parallel processing.

I received the following output:

Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/lacZ_A_S1_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/lacZ_B_S2_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/lacZ_C_S3_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/lacZ_D_S4_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/Tcf4_A_S5_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/Tcf4_B_S6_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/Tcf4_C_S7_dedub_qc_filtered.bam125 
Sample: /Volumes/day-lab$/Nathaniel/Sequencing/ATAC-seq/NJR_0019_Tcf4_CRISPRa_Striatal_ATACseq/BAM_Files/Tcf4_E_S8_dedub_qc_filtered.bam125 
Re-centering peaks...
Error in heights * sapply(called, function(x) x) : non-conformable arrays

Session info:

R version 4.2.1 (2022-06-23)
Platform: x86_64-apple-darwin17.0 (64-bit)
Running under: macOS Big Sur 11.3.1

Matrix products: default
LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

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

other attached packages:
 [1] viridis_0.6.2               viridisLite_0.4.0           pheatmap_1.0.12             ggrepel_0.9.1              
 [5] ggplot2_3.3.6               sp_1.5-0                    SeuratObject_4.1.0          Seurat_4.1.1               
 [9] reshape2_1.4.4              rtracklayer_1.56.1          DiffBind_3.6.2              SummarizedExperiment_1.26.1
[13] Biobase_2.56.0              MatrixGenerics_1.8.1        matrixStats_0.62.0          GenomicRanges_1.48.0       
[17] GenomeInfoDb_1.32.2         IRanges_2.30.0              S4Vectors_0.34.0            BiocGenerics_0.42.0        

loaded via a namespace (and not attached):
  [1] plyr_1.8.7               igraph_1.3.4             lazyeval_0.2.2           splines_4.2.1            BiocParallel_1.30.3     
  [6] listenv_0.8.0            scattermore_0.8          amap_0.8-18              digest_0.6.29            invgamma_1.1            
 [11] htmltools_0.5.3          SQUAREM_2021.1           fansi_1.0.3              magrittr_2.0.3           BSgenome_1.64.0         
 [16] tensor_1.5               cluster_2.1.3            ROCR_1.0-11              limma_3.52.2             globals_0.15.1          
 [21] Biostrings_2.64.0        systemPipeR_2.2.2        bdsmatrix_1.3-6          spatstat.sparse_2.1-1    jpeg_0.1-9              
 [26] colorspace_2.0-3         apeglm_1.18.0            dplyr_1.0.9              crayon_1.5.1             RCurl_1.98-1.8          
 [31] jsonlite_1.8.0           progressr_0.10.1         spatstat.data_2.2-0      survival_3.3-1           zoo_1.8-10              
 [36] glue_1.6.2               polyclip_1.10-0          gtable_0.3.0             zlibbioc_1.42.0          XVector_0.36.0          
 [41] leiden_0.4.2             DelayedArray_0.22.0      future.apply_1.9.0       abind_1.4-5              scales_1.2.0            
 [46] mvtnorm_1.1-3            DBI_1.1.3                spatstat.random_2.2-0    miniUI_0.1.1.1           Rcpp_1.0.9              
 [51] xtable_1.8-4             emdbook_1.3.12           reticulate_1.25          spatstat.core_2.4-4      truncnorm_1.0-8         
 [56] htmlwidgets_1.5.4        httr_1.4.3               gplots_3.1.3             RColorBrewer_1.1-3       ellipsis_0.3.2          
 [61] ica_1.0-3                pkgconfig_2.0.3          XML_3.99-0.10            uwot_0.1.11              deldir_1.0-6            
 [66] locfit_1.5-9.6           utf8_1.2.2               tidyselect_1.1.2         rlang_1.0.4              later_1.3.0             
 [71] munsell_0.5.0            tools_4.2.1              cli_3.3.0                generics_0.1.3           ggridges_0.5.3          
 [76] stringr_1.4.0            fastmap_1.1.0            goftest_1.2-3            yaml_2.3.5               fitdistrplus_1.1-8      
 [81] caTools_1.18.2           purrr_0.3.4              RANN_2.6.1               nlme_3.1-158             pbapply_1.5-0           
 [86] future_1.27.0            mime_0.12                compiler_4.2.1           plotly_4.10.0            png_0.1-7               
 [91] spatstat.utils_2.3-1     tibble_3.1.8             stringi_1.7.8            rgeos_0.5-9              lattice_0.20-45         
 [96] Matrix_1.4-1             vctrs_0.4.1              pillar_1.8.0             lifecycle_1.0.1          spatstat.geom_2.4-0     
[101] lmtest_0.9-40            RcppAnnoy_0.0.19         data.table_1.14.2        cowplot_1.1.1            bitops_1.0-7            
[106] irlba_2.3.5              httpuv_1.6.5             patchwork_1.1.1          R6_2.5.1                 BiocIO_1.6.0            
[111] latticeExtra_0.6-30      hwriter_1.3.2.1          promises_1.2.0.1         ShortRead_1.54.0         KernSmooth_2.23-20      
[116] gridExtra_2.3            parallelly_1.32.1        codetools_0.2-18         MASS_7.3-58.1            gtools_3.9.3            
[121] rjson_0.2.21             withr_2.5.0              GenomicAlignments_1.32.1 sctransform_0.3.3        Rsamtools_2.12.0        
[126] GenomeInfoDbData_1.2.8   mgcv_1.8-40              parallel_4.2.1           rpart_4.1.16             grid_4.2.1              
[131] tidyr_1.2.0              coda_0.19-4              GreyListChIP_1.28.1      ashr_2.2-54              Rtsne_0.16              
[136] mixsqp_0.3-43            bbmle_1.0.25             numDeriv_2016.8-1.1      shiny_1.7.2              interp_1.1-3            
[141] restfulr_0.0.15

One possibility I could envision is that the BAM files were somehow corrupted when I copied them over from our HPC. However, I have used these files for other analyses with no issues, so that may be unlikely, but thought I would mention it anyway. I appreciate your continued help!

ADD REPLY
0
Entering edit mode

I've checked in a fix for this (DiffBind_3.6.3).

It had to do with specifying minOverlap=0 in the initial call to dba(). Did you mean for the consensus peakset used in dba.count() to include all the (merged) peaks? If so you would need to specify minOverlap=0 in the call to dba.count() as well (or instead). Indeed, if you move the minOverlap=2 from the dba() call to the dba.count() call, your code will run fine even without the fix.

Also, you should have a very specific reason to set bUseSummarizeOverlaps=FALSE, especially if you have any paired-end sequencing data.

ADD REPLY
0
Entering edit mode

Thanks a lot Rory! I wanted to follow up for a bit more information about the potential issues encountered when setting bUseSummarizeOverlaps=FALSE with paired-end sequencing data. In my case (and others in my lab), we are analyzing ATAC-seq data (PE75). As part of our pipeline, we remove all duplicate reads using Picard. We also remove any reads that are not properly paired. Given these steps, would setting bUseSummarizeOverlaps=FALSE still be potentially problematic? And, more broadly, what are the inherent limitations when using summarizeOverlap with paired-end data? I have found some conflicting information about this in various places online, and depending on the nature of the limitations, my lab and I are curious whether we should re-analyze some of our existing data.

Thanks again for your continued help!

ADD REPLY

Login before adding your answer.

Traffic: 417 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