For a newer DiffBind version
library(DiffBind)
library(BiocParallel)
library(dplyr)
tamoxifen <- dba(sampleSheet="tamoxifen.csv")
tamoxifen <- dba.count(tamoxifen, summits=250)
compare_function <- function(group1, group2){
dba_contrast <- dba.contrast(tamoxifen,
group1 = tamoxifen$masks[[group1]],
group2 = tamoxifen$masks[[group2]],
name1 = group1,
name2 = group2)
dba_analyze <- dba.analyze(dba_contrast)
dba_report_all <- dba.report(dba_analyze,
th = 1,
precision = 0)
dba_report_all$feature_id <- paste0("peak", names(dba_report_all))
return(dba_report_all)
}
compare_1 <- compare_function(group1 = "BT474",
group2 = "MCF7")
compare_2 <- compare_function(group1 = "BT474",
group2 = "T47D")
> subset(compare_1, feature_id == "peak1006")$Conc_BT474
[1] 1.983722
> subset(compare_2, feature_id == "peak1006")$Conc_BT474
[1] 2.190726
> inner_join(data_1, data_2) %>%
+ mutate(compare = compare_2 - compare_1) %>%
+ pull(compare) %>%
+ summary()
Joining, by = "feature_id"
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.0000 0.2070 0.2070 0.2063 0.2070 0.2070
> sessionInfo()
R version 4.0.2 (2020-06-22)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19041)
Matrix products: default
locale:
[1] LC_COLLATE=Chinese (Simplified)_China.936 LC_CTYPE=Chinese (Simplified)_China.936
[3] LC_MONETARY=Chinese (Simplified)_China.936 LC_NUMERIC=C
[5] LC_TIME=Chinese (Simplified)_China.936
attached base packages:
[1] parallel stats4 stats graphics grDevices utils datasets methods base
other attached packages:
[1] dplyr_1.0.2 BiocParallel_1.24.1 DiffBind_3.0.10
[4] SummarizedExperiment_1.20.0 Biobase_2.50.0 MatrixGenerics_1.2.0
[7] matrixStats_0.57.0 GenomicRanges_1.42.0 GenomeInfoDb_1.26.2
[10] IRanges_2.24.1 S4Vectors_0.28.1 BiocGenerics_0.36.0
loaded via a namespace (and not attached):
[1] backports_1.2.0 GOstats_2.56.0 BiocFileCache_1.14.0
[4] plyr_1.8.6 GSEABase_1.52.1 splines_4.0.2
[7] ggplot2_3.3.3 amap_0.8-18 digest_0.6.27
[10] invgamma_1.1 GO.db_3.12.1 fansi_0.4.1
[13] SQUAREM_2021.1 magrittr_2.0.1 checkmate_2.0.0
[16] memoise_1.1.0 BSgenome_1.58.0 base64url_1.4
[19] limma_3.46.0 Biostrings_2.58.0 annotate_1.68.0
[22] systemPipeR_1.24.2 askpass_1.1 bdsmatrix_1.3-4
[25] prettyunits_1.1.1 jpeg_0.1-8.1 colorspace_2.0-0
[28] blob_1.2.1 rappdirs_0.3.1 apeglm_1.12.0
[31] ggrepel_0.9.0 crayon_1.3.4 RCurl_1.98-1.2
[34] jsonlite_1.7.2 graph_1.68.0 genefilter_1.72.0
[37] brew_1.0-6 survival_3.2-7 VariantAnnotation_1.36.0
[40] glue_1.4.2 gtable_0.3.0 zlibbioc_1.36.0
[43] XVector_0.30.0 DelayedArray_0.16.0 V8_3.4.0
[46] Rgraphviz_2.34.0 scales_1.1.1 pheatmap_1.0.12
[49] mvtnorm_1.1-1 DBI_1.1.0 edgeR_3.32.0
[52] Rcpp_1.0.5 xtable_1.8-4 progress_1.2.2
[55] emdbook_1.3.12 bit_4.0.4 rsvg_2.1
[58] AnnotationForge_1.32.0 truncnorm_1.0-8 httr_1.4.2
[61] gplots_3.1.1 RColorBrewer_1.1-2 ellipsis_0.3.1
[64] pkgconfig_2.0.3 XML_3.99-0.5 dbplyr_2.0.0
[67] locfit_1.5-9.4 tidyselect_1.1.0 rlang_0.4.10
[70] AnnotationDbi_1.52.0 munsell_0.5.0 tools_4.0.2
[73] cli_2.2.0 generics_0.1.0 RSQLite_2.2.2
[76] stringr_1.4.0 yaml_2.2.1 bit64_4.0.5
[79] caTools_1.18.0 purrr_0.3.4 RBGL_1.66.0
[82] xml2_1.3.2 biomaRt_2.46.0 compiler_4.0.2
[85] rstudioapi_0.13 curl_4.3 png_0.1-7
[88] geneplotter_1.68.0 tibble_3.0.4 stringi_1.5.3
[91] GenomicFeatures_1.42.1 lattice_0.20-41 Matrix_1.2-18
[94] vctrs_0.3.6 pillar_1.4.7 lifecycle_0.2.0
[97] data.table_1.13.6 bitops_1.0-6 irlba_2.3.3
[100] rtracklayer_1.49.5 R6_2.5.0 latticeExtra_0.6-29
[103] hwriter_1.3.2 ShortRead_1.48.0 KernSmooth_2.23-18
[106] MASS_7.3-53 gtools_3.8.2 assertthat_0.2.1
[109] DESeq2_1.30.0 openssl_1.4.3 Category_2.56.0
[112] rjson_0.2.20 withr_2.3.0 GenomicAlignments_1.26.0
[115] batchtools_0.9.15 Rsamtools_2.6.0 GenomeInfoDbData_1.2.4
[118] hms_0.5.3 grid_4.0.2 DOT_0.1
[121] coda_0.19-4 GreyListChIP_1.22.0 ashr_2.2-47
[124] mixsqp_0.3-43 bbmle_1.0.23.1 numDeriv_2016.8-1.1
I use
as_tibble
to convertGrange
objecet, and then use write_csv. The basic command is like belowBy the way, I found the N in
Binding Affinity: treat vs. control (N FDR < 0.050)ยท
of MA-plot is different from the N result by using filter in R or excel. I hope this phenomenon may be helpful.For example, the N in MA-plot text may be 6809 FDR < 0.05 while in my result is 6801.
dba.plotMA()
uses(abs(Fold) >= fold) & (FDR <= th)
.For the concentration issue, if you could send me a link where I can access your
dba_diff
object, I can see if I can reproduce it. I'll also need to see the output from callingsessionInfo()
so I know which software versions and platform you are running.I am sorry I did not post link and
seesionInfo
Here is the link, I help it will be helpful
And the session
Hi, Dr Stark
I just find you fixed the link, and I want to use the test data to reproduce my issue