get_ssRleCov {InPAS} | R Documentation |
Get RLe coverage from a bedgraph file for a sample
get_ssRleCov( bedgraph, tag, genome, sqlite_db, outdir, BATCH_SIZE = 10L, removeScaffolds = FALSE, BPPARAM = NULL )
bedgraph |
A path to a bedGraph file |
tag |
A character(1) vector, a name tag used to label the bedgraph file. It must match the tag specified in the metadata file used to setup the SQLite database |
genome |
an object BSgenome::BSgenome. To make things easy, we
suggest users creating a BSgenome::BSgenome instance from the
reference genome used for read alignment. For details, see the
documentation of |
sqlite_db |
A path to the SQLite database for InPAS, i.e. the output of
|
outdir |
A character(1) vector, a path with write permission for storing the coverage data. If it doesn't exist, it will be created. |
BATCH_SIZE |
A integer(1) vector, indicating the number of parallel jobs run at the same time per batch. Default, 10. You may adjust this number based based on the available computing resource: CPUs and RAM. For BATCH_SIZE of 10, 15-20G RAM is needed. This parameter affects the time for converting coverage from bedgraph to Rle. |
removeScaffolds |
A logical(1) vector, whether the scaffolds should be removed from the genome If you use a TxDb containing alternative scaffolds, you'd better to remove the scaffolds. |
BPPARAM |
an optional BiocParallel::BiocParallelParam instance determining the parallel back-end to be used during evaluation, or a list of BiocParallelParam instances, to be applied in sequence for nested calls to bplapply. It can be set to NULL or bpparam() |
A list of lists containing read coverage as Rle instances of S4Vectors::Rle representing read coverage for each chromosome of a given sample, as described below.
the sample tag
coverage as Rle instance for chr1
coverage as Rle instance for chr2
coverage as Rle instance for chrN
Jianhong Ou, Haibo Liu
if (interactive()) { library(BSgenome.Mmusculus.UCSC.mm10) genome <- BSgenome.Mmusculus.UCSC.mm10 bedgraphs <- system.file("extdata",c("Baf3.extract.bedgraph", "UM15.extract.bedgraph"), package = "InPAS") tags <- c("Baf3", "UM15") metadata <- data.frame(tag = tags, condition = c("Baf3", "UM15"), bedgraph_file = bedgraphs) outdir = tempdir() write.table(metadata, file =file.path(outdir, "metadata.txt"), sep = "\t", quote = FALSE, row.names = FALSE) sqlite_db <- setup_sqlitedb(metadata = file.path(outdir, "metadata.txt"), outdir) coverage <- get_ssRleCov(bedgraph = bedgraphs[1], tag = tags[1], genome = genome, sqlite_db = sqlite_db, outdir = outdir, removeScaffolds = TRUE, BPPARAM = NULL) # check read coverage depth db_connect <- dbConnect(drv = RSQLite::SQLite(), dbname = sqlite_db) dbReadTable(db_connect, "metadata") dbDisconnect(db_connect) }