I am running
derfinder to ultimately find DER in a time course but despite the very good documentation I have some doubts about how set the cutoff to find expressed regions. In the 2017 NAR paper, it is stated that an average across all regions and samples is used as a cut-off (Finding ERs):
Here, we focus on a new approach based on the bump-hunting methodology for region level genomic analysis (33) (Figure 2). This approach first calculates ERs across the set of observed samples. For each base, the average, potentially library size-adjusted, coverage is calculated across all samples in the data set. This generates a vector of (normalized) mean level expression measurements across the genome. Then an average-coverage cutoff is applied to this mean coverage vector to identify bases that show minimum levels of expression. An expressed region is any contiguous set of bases that has expression above the mean expression cutoff.
However, in the manual I could not find any indication of this being used, and all examples use a hard cutoff for base coverage. I also looked into help page of
filterData but I am still in the dark.
- Is this filtering option using calculated mean coverage implemented? I could of course calculate that myself and feed it to
regionMatrixbut that could be re-inventing the wheel.
- If it turns out I am misunderstanding the above paper paragraph, any tips on selecting a reasonable cut-off? I would usually create a distribution of random bases (sampling n times) and then use the median as a threshold.
Unrelated to the main points, the data are normalized bigwigs, generated from paired-end sequencing. Which read-length should bed used in
regionMatrix? Actual read length or something else? derfinder on bacteria PE strand specific RNA-seq is somewhat related but I am not interested in finding expressed regions for each individual strand (at least not at the moment).
> sessionInfo() R version 3.4.1 (2017-06-30) Platform: x86_64-pc-linux-gnu (64-bit) Running under: Debian GNU/Linux 8 (jessie) Matrix products: default BLAS: /usr/lib/atlas-base/libf77blas.so.3.0 LAPACK: /usr/lib/atlas-base/atlas/liblapack.so.3.0 locale:  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  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:  parallel stats4 stats graphics grDevices utils datasets  methods base other attached packages:  DESeq2_1.18.1 SummarizedExperiment_1.8.0  DelayedArray_0.4.1 matrixStats_0.52.2  Biobase_2.38.0 GenomicRanges_1.30.0  GenomeInfoDb_1.14.0 IRanges_2.12.0  S4Vectors_0.16.0 BiocGenerics_0.24.0  derfinder_1.12.0 fortunes_1.5-4 loaded via a namespace (and not attached):  RMySQL_0.10.13 bit64_0.9-7 splines_3.4.1  foreach_1.4.3 GenomicFiles_1.14.0 Formula_1.2-2  bumphunter_1.20.0 assertthat_0.2.0 latticeExtra_0.6-28  doRNG_1.6.6 blob_1.1.0 BSgenome_1.46.0  GenomeInfoDbData_0.99.1 Rsamtools_1.30.0 progress_1.1.2  RSQLite_2.0 backports_1.1.1 lattice_0.20-35  digest_0.6.12 RColorBrewer_1.1-2 XVector_0.18.0  checkmate_1.8.5 qvalue_2.10.0 colorspace_1.3-2  htmltools_0.3.6 Matrix_1.2-12 plyr_1.8.4  devtools_1.13.4 XML_3.98-1.9 biomaRt_2.34.0  genefilter_1.60.0 zlibbioc_1.24.0 xtable_1.8-2  scales_0.5.0 BiocParallel_1.12.0 annotate_1.56.1  htmlTable_1.9 tibble_1.3.4 pkgmaker_0.22  ggplot2_2.2.1 withr_2.1.0 GenomicFeatures_1.30.0  nnet_7.3-12 lazyeval_0.2.1 survival_2.41-3  magrittr_1.5 memoise_1.1.0 foreign_0.8-69  registry_0.3 tools_3.4.1 data.table_1.10.4-3  prettyunits_1.0.2 stringr_1.2.0 locfit_1.5-9.1  munsell_0.4.3 rngtools_1.2.4 cluster_2.0.6  AnnotationDbi_1.40.0 Biostrings_2.46.0 compiler_3.4.1  rlang_0.1.4 grid_3.4.1 RCurl_1.95-4.8  rstudioapi_0.7 iterators_1.0.8 VariantAnnotation_1.24.1  htmlwidgets_0.9 bitops_1.0-6 base64enc_0.1-3  derfinderHelper_1.12.0 gtable_0.2.0 codetools_0.2-15  DBI_0.7 reshape2_1.4.2 R6_2.2.2  GenomicAlignments_1.14.0 gridExtra_2.3 knitr_1.17  rtracklayer_1.38.0 bit_1.1-12 Hmisc_4.0-3  stringi_1.1.6 Rcpp_0.12.13 geneplotter_1.56.0  rpart_4.1-11 acepack_1.4.1