I have run ATACSeq on 5 sample: 3 case two control. I called the peaks by genrich by merging case and control sampes. I am having two narropeak files. Now, I am trying to use diffbind for identifying differtial peaks. But I am have the following error whole running dba.count. I am not sure why I am having this. Any suggestions?
Thank you
dbObj <- dba.count(dbObj, bUseSummarizeOverlaps = TRUE) Computing summits... Re-centering peaks... Reads will be counted as Paired-end.
Error in DESeqDataSet(se, design = design, ignoreRank) : all samples have 0 counts for all genes. check the counting script.
sampleSheet
SampleID,Tissue,Factor,Condition,Replicate,Treatment,bamReads,ControlID,bamControl,Peaks,PeakCaller
C1,CK,CE,Case,1,Full-Media,C1.bam,CE-1,C1.bam,sample.narroPeak,narrow
C2,CK,CE,Case,2,Full-Media,C2.bam,CE-2,C2.bam,sample.narroPeak,narrow
C3,CK,CE,Case,3,Full-Media,C3.bam,CE-3,C3.bam,sample.narroPeak,narrow
HT_1,HET,HK,Control,1,Full-Media,HT_1.bam,HK-1,HT_1.bam,control.narroPeak,narrow
HT_2,HET,HK,Control,2,Full-Media,HT_2.bam,HK-2,HT_2.bam,control.narroPeak,narrow`
sessionInfo( )
``` R version 4.0.3 (2020-10-10) Platform: x86_64-apple-darwin17.0 (64-bit) Running under: macOS Mojave 10.14
Matrix products: default BLAS: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRblas.dylib LAPACK: /Library/Frameworks/R.framework/Versions/4.0/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] parallel stats4 stats graphics grDevices utils datasets [8] methods base
other attached packages:
[1] forcats_0.5.0 stringr_1.4.0
[3] dplyr_1.0.2 purrr_0.3.4
[5] readr_1.4.0 tidyr_1.1.2
[7] tibble_3.0.4 ggplot2_3.3.2
[9] tidyverse_1.3.0 DiffBind_3.0.6
[11] SummarizedExperiment_1.20.0 Biobase_2.50.0
[13] MatrixGenerics_1.2.0 matrixStats_0.57.0
[15] GenomicRanges_1.42.0 GenomeInfoDb_1.26.0
[17] IRanges_2.24.0 S4Vectors_0.28.0
[19] BiocGenerics_0.36.0
loaded via a namespace (and not attached):
[1] readxl_1.3.1 backports_1.2.0 GOstats_2.56.0
[4] BiocFileCache_1.14.0 plyr_1.8.6 GSEABase_1.52.0
[7] splines_4.0.3 BiocParallel_1.24.1 amap_0.8-18
[10] digest_0.6.27 invgamma_1.1 GO.db_3.12.1
[13] fansi_0.4.1 SQUAREM_2020.5 magrittr_2.0.1
[16] checkmate_2.0.0 memoise_1.1.0 BSgenome_1.58.0
[19] base64url_1.4 limma_3.46.0 Biostrings_2.58.0
[22] annotate_1.68.0 modelr_0.1.8 systemPipeR_1.24.2
[25] askpass_1.1 bdsmatrix_1.3-4 prettyunits_1.1.1
[28] jpeg_0.1-8.1 colorspace_2.0-0 rvest_0.3.6
[31] blob_1.2.1 rappdirs_0.3.1 apeglm_1.12.0
[34] ggrepel_0.8.2 haven_2.3.1 crayon_1.3.4
[37] RCurl_1.98-1.2 jsonlite_1.7.1 graph_1.68.0
[40] genefilter_1.72.0 brew_1.0-6 survival_3.2-7
[43] VariantAnnotation_1.36.0 glue_1.4.2 gtable_0.3.0
[46] zlibbioc_1.36.0 XVector_0.30.0 DelayedArray_0.16.0
[49] V8_3.4.0 Rgraphviz_2.34.0 scales_1.1.1
[52] pheatmap_1.0.12 mvtnorm_1.1-1 DBI_1.1.0
[55] edgeR_3.32.0 Rcpp_1.0.5 xtable_1.8-4
[58] progress_1.2.2 emdbook_1.3.12 bit_4.0.4
[61] rsvg_2.1 AnnotationForge_1.32.0 truncnorm_1.0-8
[64] httr_1.4.2 gplots_3.1.0 RColorBrewer_1.1-2
[67] ellipsis_0.3.1 pkgconfig_2.0.3 XML_3.99-0.5
[70] dbplyr_2.0.0 locfit_1.5-9.4 tidyselect_1.1.0
[73] rlang_0.4.8 AnnotationDbi_1.52.0 cellranger_1.1.0
[76] munsell_0.5.0 tools_4.0.3 cli_2.1.0
[79] generics_0.1.0 RSQLite_2.2.1 broom_0.7.2
[82] yaml_2.2.1 fs_1.5.0 bit64_4.0.5
[85] caTools_1.18.0 RBGL_1.66.0 xml2_1.3.2
[88] biomaRt_2.46.0 rstudioapi_0.13 compiler_4.0.3
[91] curl_4.3 png_0.1-7 reprex_0.3.0
[94] geneplotter_1.68.0 stringi_1.5.3 ps_1.4.0
[97] GenomicFeatures_1.42.1 lattice_0.20-41 Matrix_1.2-18
[100] vctrs_0.3.5 pillar_1.4.6 lifecycle_0.2.0
[103] data.table_1.13.2 bitops_1.0-6 irlba_2.3.3
[106] rtracklayer_1.50.0 R6_2.5.0 latticeExtra_0.6-29
[109] hwriter_1.3.2 ShortRead_1.48.0 KernSmooth_2.23-18
[112] MASS_7.3-53 gtools_3.8.2 assertthat_0.2.1
[115] DESeq2_1.30.0 openssl_1.4.3 Category_2.56.0
[118] rjson_0.2.20 withr_2.3.0 GenomicAlignments_1.26.0
[121] batchtools_0.9.14 Rsamtools_2.6.0 GenomeInfoDbData_1.2.4
[124] hms_0.5.3 grid_4.0.3 DOT_0.1
[127] coda_0.19-4 GreyListChIP_1.22.0 ashr_2.2-47
[130] mixsqp_0.3-43 bbmle_1.0.23.1 lubridate_1.7.9.2
[133] numDeriv_2016.8-1.1
NB: This issue was resolved in email.
Hi, I am having the same issue. Was wondering if you could explain how it was resolved. Thank you!
In that case, the user was specifying the same ChIP/ATAC read files as both the signal and control reads, and not using a greylist. By default, in this case, the control reads are subtracted form the signal reads, resulting in all zero counts.
In other cases, users have gotten all zero counts when for some reason the chromosome names listed in the peak files were different than those used in the bam files.
If neither of these fits your case, you can open a new issue here with more details as to your script and I'll see if I can help.
I am Xiaolong Ma, a bioinformatics analyst from Medical College of Wisconsin, the U.S.
Recently, I have bumped against a tough issue when running the Diffbind package, it might be related to the problem above, which was that I could pass the dba.count function and get the results, but the result did not have FRiP column. I have clipped my code and paste below. Could you please advise me how to deal with this problem? Thank you!
Counting reads:
db0323yw<- dba.count(db0323yw, bUseSummarizeOverlaps=TRUE, bParallel=FALSE)
plot(db0323yw);db0323yw
1st row(overlap in at least two of the samples)
2nd row(total number of unique peaks after merging overlapping ones)
FRiP=Fraction of Reads in Peaks
info <- dba.show(db0323yw) libsizes <- cbind(LibReads=info$Reads, FRiP=info$FRiP,PeakReads=round(info$Reads * info$FRiP)) rownames(libsizes) <- info$ID write.table(libsizes,"libsizes.txt",quote=F)
saveRDS(db0323yw, file="db0323yw_2.rds")
load("db0323yw.rds")
Establishing a contrast: groupby according to factor, 3 groups
db0323yw<- dba.contrast(db0323yw, categories=DBA_FACTOR, minMembers = 3)
BTW, I had this error when running contrast function.