diff --git a/bin/nucDyn.R b/bin/nucDyn.R index f296314..7751635 100644 --- a/bin/nucDyn.R +++ b/bin/nucDyn.R @@ -13,6 +13,7 @@ suppressPackageStartupMessages(library(NucDyn)) suppressPackageStartupMessages(library(plyr)) suppressPackageStartupMessages(library(parallel)) + where <- function () { spath <- parent.frame(2)$ofile @@ -126,7 +127,6 @@ hs <- findHotspots(dyn=dyn, nuc=nuc, mc.cores=params$cores , ## Calculate vector of -log10(p-value)s ###################################### message("calculate p-val") - cov1 <- lapply(coverage(as(rs[[1]], "GRanges")), as.vector) cov2 <- lapply(coverage(as(rs[[2]], "GRanges")), as.vector) @@ -134,7 +134,6 @@ chrs <- intersect(names(cov1), names(cov2)) pvals <- lapply(chrs, function (x) -log10(findPVals(cov1[[x]], cov2[[x]]))) names(pvals) <- chrs - ## Store the Result ########################################################### names(hs)[names(hs) == "type"] <- "class" @@ -148,10 +147,15 @@ writeGff(df2gff(hs, params$outputGff) message("saving bigWig output") -writeBigWig(lapply(pvals, splitAtZeros), - params$outputBigWig, - params$genome) + +writeBigWig_updated(pvals, + params$outputBigWig, + params$genome) + +# writeBigWig(lapply(pvals, splitAtZeros), +# params$outputBigWig, +# params$genome) message("done") -############################################################################### +############################################################################### \ No newline at end of file diff --git a/bin/periodicity.R b/bin/periodicity.R index eb0eb09..52925cc 100644 --- a/bin/periodicity.R +++ b/bin/periodicity.R @@ -114,7 +114,10 @@ gff <- df2gff(genes.nucs, writeGff(gff, params$gffOutput) message("writing bigWig output") -splited <- lapply(covPredAll, splitAtZeros) -writeBigWig(splited, params$bwOutput, params$chrom_sizes) +# splited <- lapply(covPredAll, splitAtZeros) +# writeBigWig(splited, params$bwOutput, params$chrom_sizes) +writeBigWig_updated(covPredAll, + params$bwOutput, + params$chrom_sizes) ############################################################################## diff --git a/sourced/wig_funs.R b/sourced/wig_funs.R index 9e4e1a7..b6890b1 100644 --- a/sourced/wig_funs.R +++ b/sourced/wig_funs.R @@ -3,6 +3,10 @@ ## Imports #################################################################### suppressPackageStartupMessages(library(IRanges)) +suppressPackageStartupMessages(library(plyranges)) +suppressPackageStartupMessages(library(purrr)) +suppressPackageStartupMessages(library(dplyr)) +suppressPackageStartupMessages(library(readr)) source(paste(SOURCE.DIR, "get_genes.R", sep="/")) ## Binary paths ############################################################### @@ -118,5 +122,39 @@ writeBigWig <- function (x, outf, chrom.sizes.f) file.remove(wigf) invisible() } +# This function takes a list of vectors (x) representing signal scores across genomic regions and writes the data to a bigWig file format +# Additionally it uses a file path to write the resulting bigWig file (outf), +# and a file path to a file with chromosome sizes in tab seperated form (chrom.sizes.f) +writeBigWig_updated <- function (x, outf, chrom.sizes.f) { + + # Map over the input vectors to create a tibble for each vector with four columns: + # score, start position, width, and strand (* indicates unstranded). + # Combine all tibbles into a single data frame using map_df(). + # Use .id = "seqnames" to create a new column seqnames, which contains the names of the vectors in x. + gr <- purrr::map_df(x, + function(y) + {dplyr::tibble(score = y, + start = 1:length(y), + width = 1, + strand = "*")}, + .id = "seqnames") %>% + as_granges() + + # Read in a file containing chromosome sizes into a data frame. + # The file is assumed to have two columns separated by tabs: chromosome name and size. + chrSizes <- readr::read_delim(chrom.sizes.f, + col_names = c("seqnames","size"), + delim = "\t", + show_col_types = FALSE) + + # Filter the chromosome sizes data frame to include only the chromosomes that are present in the granges object. + chrSizes <- filter(chrSizes, seqnames %in% seqnames(gr)) + + # Set the sequence lengths of the granges object to the chromosome sizes in chrSizes. + seqlengths(gr) <- chrSizes$size + + # Write the granges object to a bigWig file at the specified output path. + write_bigwig(gr, outf) +} ###############################################################################