diff --git a/DESCRIPTION b/DESCRIPTION index 87e55b16..46f19b5a 100755 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Type: Package Package: misha Title: Toolkit for Analysis of Genomic Data -Version: 5.11.11 +Version: 5.11.12 Authors@R: c( person("Misha", "Hoichman", , "misha@hoichman.com", role = "aut"), person("Aviezer", "Lifshitz", , "aviezer.lifshitz@weizmann.ac.il", role = c("aut", "cre")), diff --git a/NEWS.md b/NEWS.md index 496ded4e..d84e1c5a 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,8 @@ +# misha 5.11.12 + +* `gtrack.import()` no longer fails on non-Linux platforms: the bundled `bigWigToWig` is a Linux x86-64 binary, so macOS/Windows now use one found on `PATH` (or `options(misha.bigWigToWig = "...")`), with an actionable error if none is installed. +* `gtrack.create_sparse()` can take the values from a `value` column of `intervals` when `values` is omitted, which rules out a values/intervals order mismatch. Documented that `gintervals()` returns its rows sorted, so a separately held value vector must be reordered too. + # misha 5.11.11 * **Internal consistency:** `N`/`*` are now scored the same on both strands in the PWM energy routines (reverse used `log(0.25)`, forward used the column average; both now use the average). No behavior change for genomic scoring, where `N`-windows are masked to `-Inf`; keeps `DnaPSSM` in sync with the `prego` package. diff --git a/R/intervals-core.R b/R/intervals-core.R index 2534f32c..46a7203c 100644 --- a/R/intervals-core.R +++ b/R/intervals-core.R @@ -105,7 +105,11 @@ #' @param ends an array of end coordinates. If '-1' chromosome size is assumed. #' @param strands 'NULL', a numeric vector of '-1', '0' or '1' values, or a #' character/factor vector with values "+", "-", ".", "*" or "" -#' @return A data frame representing the intervals. +#' @return A data frame representing the intervals, sorted in the canonical +#' chromosome order (which is not the order of the arguments). Beware of keeping +#' per-interval data in a separate parallel vector across a \code{gintervals()} +#' call: the rows are reordered and the vector is not. Put such data in a column +#' of the resulting data frame instead. #' @seealso \code{\link{gintervals.2d}}, \code{\link{gintervals.force_range}} #' @keywords ~intervals #' @examples diff --git a/R/track-create.R b/R/track-create.R index b38326f2..49a43a50 100644 --- a/R/track-create.R +++ b/R/track-create.R @@ -215,10 +215,20 @@ gtrack.create_pwm_energy <- function(track = NULL, description = NULL, pssmset = #' last connected database. Use \code{\link{gdir.cd}} with an absolute path to #' change where new tracks are created. #' +#' \code{values} is matched to \code{intervals} row by row, in the order the +#' intervals are passed; \code{intervals} need not be sorted. Note however that +#' \code{\link{gintervals}} returns its result sorted in the canonical +#' chromosome order, so building \code{intervals} with \code{gintervals()} while +#' keeping \code{values} in the original order will bind values to the wrong +#' intervals. Keep the values in a \code{value} column of the intervals data +#' frame and omit \code{values} to make such a mismatch impossible. +#' #' @param track track name #' @param description a character string description #' @param intervals a set of one-dimensional intervals -#' @param values an array of numeric values - one for each interval +#' @param values an array of numeric values - one for each interval, in the same +#' order as the rows of \code{intervals}. If \code{NULL}, the \code{value} +#' column of \code{intervals} is used. #' @return None. #' @seealso \code{\link{gtrack.create}}, \code{\link{gtrack.2d.create}}, #' \code{\link{gtrack.smooth}}, \code{\link{gtrack.modify}}, @@ -241,13 +251,20 @@ gtrack.create_pwm_energy <- function(track = NULL, description = NULL, pssmset = #' #' @export gtrack.create_sparse gtrack.create_sparse <- function(track = NULL, description = NULL, intervals = NULL, values = NULL) { - if (is.null(substitute(track)) || is.null(description) || is.null(intervals) || is.null(values)) { + if (is.null(substitute(track)) || is.null(description) || is.null(intervals)) { stop("Usage: gtrack.create_sparse(track, description, intervals, values)", call. = FALSE) } .gcheckroot() intervals <- rescue_ALLGENOME(intervals, as.character(substitute(intervals))) + if (is.null(values)) { + if (!is.data.frame(intervals) || !("value" %in% colnames(intervals))) { + stop("values is missing and intervals has no \"value\" column", call. = FALSE) + } + values <- intervals$value + } + trackstr <- do.call(.gexpr2str, list(substitute(track)), envir = parent.frame()) .gconfirmtrackcreate(trackstr) diff --git a/R/track-import.R b/R/track-import.R index a3fe6c9b..78f06f03 100644 --- a/R/track-import.R +++ b/R/track-import.R @@ -1,9 +1,29 @@ # Track import functions and helpers -get_bigWigToWig_bin <- function() { - dir <- tempdir() - utils::untar(system.file("bigWigToWig.tar.gz", package = "misha"), exdir = dir) - return(file.path(dir, "bigWigToWig")) +get_bigWigToWig_bin <- function(sysname = Sys.info()[["sysname"]]) { + bin <- getOption("misha.bigWigToWig", "") + if (nzchar(bin)) { + return(bin) + } + + # ponytail: only a Linux x86-64 binary is bundled; elsewhere fall back to PATH + if (sysname == "Linux") { + dir <- tempdir() + utils::untar(system.file("bigWigToWig.tar.gz", package = "misha"), exdir = dir) + return(file.path(dir, "bigWigToWig")) + } + + bin <- Sys.which("bigWigToWig") + if (!nzchar(bin)) { + stop(sprintf(paste0( + "The bigWigToWig executable bundled with misha is a Linux x86-64 binary and cannot run on %s.\n", + "Install UCSC's bigWigToWig for your platform, e.g.:\n", + " conda install -c bioconda ucsc-bigwigtowig\n", + " or download it from http://hgdownload.soe.ucsc.edu/admin/exe/\n", + "Then put it on PATH, or set options(misha.bigWigToWig = \"/path/to/bigWigToWig\")." + ), sysname), call. = FALSE) + } + bin } diff --git a/man/gintervals.Rd b/man/gintervals.Rd index ca811887..d9140881 100644 --- a/man/gintervals.Rd +++ b/man/gintervals.Rd @@ -18,7 +18,11 @@ prefixes or an array of integers (like: '1' for "chr1")} character/factor vector with values "+", "-", ".", "*" or ""} } \value{ -A data frame representing the intervals. +A data frame representing the intervals, sorted in the canonical +chromosome order (which is not the order of the arguments). Beware of keeping +per-interval data in a separate parallel vector across a \code{gintervals()} +call: the rows are reordered and the vector is not. Put such data in a column +of the resulting data frame instead. } \description{ Creates a set of 1D intervals. diff --git a/man/gtrack.create_sparse.Rd b/man/gtrack.create_sparse.Rd index 1c54277d..7503623b 100644 --- a/man/gtrack.create_sparse.Rd +++ b/man/gtrack.create_sparse.Rd @@ -18,7 +18,9 @@ gtrack.create_sparse( \item{intervals}{a set of one-dimensional intervals} -\item{values}{an array of numeric values - one for each interval} +\item{values}{an array of numeric values - one for each interval, in the same +order as the rows of \code{intervals}. If \code{NULL}, the \code{value} +column of \code{intervals} is used.} } \value{ None. @@ -34,6 +36,14 @@ When multiple databases are connected via \code{\link{gsetroot}}, the track is created in the current working directory (.misha$GWD), which defaults to the last connected database. Use \code{\link{gdir.cd}} with an absolute path to change where new tracks are created. + +\code{values} is matched to \code{intervals} row by row, in the order the +intervals are passed; \code{intervals} need not be sorted. Note however that +\code{\link{gintervals}} returns its result sorted in the canonical +chromosome order, so building \code{intervals} with \code{gintervals()} while +keeping \code{values} in the original order will bind values to the wrong +intervals. Keep the values in a \code{value} column of the intervals data +frame and omit \code{values} to make such a mismatch impossible. } \examples{ \dontshow{ diff --git a/tests/testthat/test-bigwigtowig-resolution.R b/tests/testthat/test-bigwigtowig-resolution.R new file mode 100644 index 00000000..d83ca83d --- /dev/null +++ b/tests/testthat/test-bigwigtowig-resolution.R @@ -0,0 +1,36 @@ +# The bundled bigWigToWig is a Linux x86-64 binary, so bigWig import cannot work +# on other platforms without an externally installed converter. There is no +# bigWig import test, which is why this went unnoticed on the macOS CI runner. + +test_that("options(misha.bigWigToWig) wins over everything else", { + withr::local_options(misha.bigWigToWig = "/some/where/bigWigToWig") + expect_equal(get_bigWigToWig_bin(), "/some/where/bigWigToWig") + # honored on every platform, including the bundled-binary one + expect_equal(get_bigWigToWig_bin("Linux"), "/some/where/bigWigToWig") +}) + +test_that("the bundled binary is used on Linux", { + skip_if_not(Sys.info()[["sysname"]] == "Linux") + withr::local_options(misha.bigWigToWig = NULL) + expect_true(file.exists(get_bigWigToWig_bin("Linux"))) +}) + +test_that("a converter on PATH is used off-Linux", { + withr::local_options(misha.bigWigToWig = NULL) + dir <- withr::local_tempdir() + fake <- file.path(dir, "bigWigToWig") + writeLines("#!/bin/sh\nexit 0", fake) + Sys.chmod(fake, "0755") + withr::local_envvar(PATH = dir) + + expect_equal(normalizePath(get_bigWigToWig_bin("Darwin")), normalizePath(fake)) +}) + +test_that("off-Linux with no converter errors with install instructions", { + withr::local_options(misha.bigWigToWig = NULL) + withr::local_envvar(PATH = withr::local_tempdir()) # empty -> Sys.which finds nothing + + expect_error(get_bigWigToWig_bin("Darwin"), "cannot run on Darwin") + expect_error(get_bigWigToWig_bin("Darwin"), "ucsc-bigwigtowig") + expect_error(get_bigWigToWig_bin("Darwin"), "misha\\.bigWigToWig") +}) diff --git a/tests/testthat/test-gtrack.create_sparse-value-alignment.R b/tests/testthat/test-gtrack.create_sparse-value-alignment.R new file mode 100644 index 00000000..29075fca --- /dev/null +++ b/tests/testthat/test-gtrack.create_sparse-value-alignment.R @@ -0,0 +1,54 @@ +create_isolated_test_db() + +# Reported by a user who built intervals with gintervals() (which sorts) while +# keeping values in the pre-sort order. create_sparse itself binds values by the +# original row index, so unsorted input is fine. +test_that("create_sparse binds values by original row order, not sorted order", { + tmptrack <- paste0("test.tmptrack_", sample(1:1e9, 1)) + gtrack.rm(tmptrack, force = TRUE) + withr::defer(gtrack.rm(tmptrack, force = TRUE)) + + # deliberately out of canonical chrom order + intervals <- data.frame( + chrom = c("chr2", "chr1", "chr2", "chr1"), + start = c(100, 200, 300, 400), + end = c(101, 201, 301, 401), + stringsAsFactors = FALSE + ) + values <- c(2.1, 1.1, 2.3, 1.4) + + gtrack.create_sparse(tmptrack, "Test", intervals, values) + + r <- gextract(tmptrack, gintervals.all(), colnames = "value") + expected <- intervals[order(intervals$chrom, intervals$start), ] + expect_equal(as.character(r$chrom), expected$chrom) + expect_equal(r$start, expected$start) + expect_equal(r$value, values[order(intervals$chrom, intervals$start)], tolerance = 1e-6) +}) + +test_that("create_sparse takes values from the value column when values is omitted", { + tmptrack <- paste0("test.tmptrack_", sample(1:1e9, 1)) + gtrack.rm(tmptrack, force = TRUE) + withr::defer(gtrack.rm(tmptrack, force = TRUE)) + + intervals <- data.frame( + chrom = c("chr2", "chr1", "chr2"), + start = c(100, 200, 300), + end = c(101, 201, 301), + value = c(2.1, 1.1, 2.3), + stringsAsFactors = FALSE + ) + + gtrack.create_sparse(tmptrack, "Test", intervals) + + r <- gextract(tmptrack, gintervals.all(), colnames = "value") + expect_equal(r$value, c(1.1, 2.1, 2.3), tolerance = 1e-6) +}) + +test_that("create_sparse errors when values is omitted and there is no value column", { + intervals <- gintervals(1, 100, 200) + expect_error( + gtrack.create_sparse("test.no_value_col", "Test", intervals), + "no \"value\" column" + ) +}) diff --git a/tests/testthat/test-gvtrack-clear.R b/tests/testthat/test-gvtrack-clear.R index e2e9bf0a..dea3c27c 100644 --- a/tests/testthat/test-gvtrack-clear.R +++ b/tests/testthat/test-gvtrack-clear.R @@ -1,3 +1,5 @@ +create_isolated_test_db() + test_that("gvtrack.clear removes all virtual tracks for the current working directory", { gvtrack.create("vt_clear_1", "test.fixedbin", "max") gvtrack.create("vt_clear_2", "test.fixedbin", "avg")