From 3cd633d2c796bc8dd7c6d5d9735f4aacce8abf7f Mon Sep 17 00:00:00 2001 From: Ivo Kwee Date: Sat, 1 Aug 2026 20:16:29 +0200 Subject: [PATCH 1/7] add new outlier method --- R/pgx-outlier.R | 58 ++++++++++++++++++++++++++++++++++++------------- 1 file changed, 43 insertions(+), 15 deletions(-) diff --git a/R/pgx-outlier.R b/R/pgx-outlier.R index a126ae09..95819d3b 100644 --- a/R/pgx-outlier.R +++ b/R/pgx-outlier.R @@ -4,25 +4,31 @@ ## #' @export -detectOutlierSamples <- function(X, plot = TRUE, par = NULL) { +detectOutlierSamples <- function(X, + methods = c("z.correlation", "z.distance", + "z.features", "z.isoforest")[1:3], + col=col, plot = TRUE, par = TRUE) { + + if(is.null(methods)) { + methods = c("z.correlation", "z.distance", "z.features", "z.isoforest") + } + ## correlation and distance X <- head(X[order(-matrixStats::rowSds(X, na.rm = TRUE)), ], 1000) X <- X - median(X, na.rm = TRUE) - corX <- HiClimR::fastCor(X, optBLAS = TRUE) + #corX <- HiClimR::fastCor(X, optBLAS = TRUE) + corX <- cor(X, use="pairwise.complete") distX <- as.matrix(dist(t(X))) + z1=z2=z3=z4=NULL + ## z-score based on correlation - ## cor.min - ## cor.max <- apply(abs(corX), 1, max, na.rm = TRUE) cor.median <- apply(abs(corX), 1, median, na.rm = TRUE) - ## cor.q10 <- apply(abs(corX), 1, quantile, probs = 0.1, na.rm = TRUE) x1 <- (cor.median - mean(cor.median, na.rm = TRUE)) z1 <- abs(x1 - median(x1, na.rm = TRUE)) / mad(x1, na.rm = TRUE) ## z-score based on euclidean distance dist.max <- apply(distX, 1, max, na.rm = TRUE) - ## dist.median <- apply(distX, 1, median, na.rm = TRUE) - ## dist.q90 <- apply(distX, 1, quantile, probs = 0.9, na.rm = TRUE) dist.q10 <- apply(distX, 1, quantile, probs = 0.1, na.rm = TRUE) dist.r <- dist.q10 / dist.max z2 <- abs(dist.r - median(dist.r, na.rm = TRUE)) / mad(dist.r, na.rm = TRUE) @@ -32,31 +38,53 @@ detectOutlierSamples <- function(X, plot = TRUE, par = NULL) { xz <- colMeans(xz, na.rm = TRUE) z3 <- abs(xz - median(xz, na.rm = TRUE)) / mad(xz, na.rm = TRUE) - Z <- cbind(z1, z2, z3) - colnames(Z) <- c("z.correlation", "z.distance", "z.features") + ## isoforest z-score + z4 <- outlier.isoforest_zscore(X, scale=TRUE, ndim=2, ntrees=10000) + + Z <- cbind(z1, z2, z3, z4) + colnames(Z) <- c("z.correlation", "z.distance", "z.features", "z.isoforest") + Z <- Z[,which(colnames(Z) %in% methods)] + zz <- rowMeans(Z, na.rm = TRUE) z0 <- 0.1 * mean(Z, na.rm = TRUE) zz2 <- exp(rowMeans(log(Z + z0), na.rm = TRUE)) - z0 res <- list(z.outlier = zz, z.outlier2 = zz2, Z = Z) - if (plot) plotOutlierScores(res, par = par) + if (plot) plotOutlierScores(res, par = par, col=col) return(res) } #' @export -plotOutlierScores <- function(res.outliers, z.threshold = c(3, 6, 9), par = TRUE) { +plotOutlierScores <- function(res.outliers, z.threshold = c(3, 6, 9), + col = "grey70", par = TRUE) { if (par) par(mfrow = c(2, 3), mar = c(8, 4, 2, 2)) Z <- res.outliers$Z zz <- res.outliers$z.outlier zz2 <- res.outliers$z.outlier2 barplot2 <- function(x, ...) { - barplot(x, ylim = c(0, max(10, max(Z))), ylab = "z-score", ...) + barplot(x, col = col, las = 3, + ylim = c(0, max(10, max(Z))), + ylab = "z-score", ...) abline(h = z.threshold, lty = 3, col = "red") } - barplot2(zz, main = "z.outlier (mean)", las = 3) - barplot2(zz2, main = "z.outlier (geom.mean)", las = 3) + barplot2(zz, main = "z.outlier (mean)") + barplot2(zz2, main = "z.outlier (geom.mean)") for (i in 1:ncol(Z)) { z1 <- Z[, i] - barplot2(z1, main = colnames(Z)[i], las = 3) + barplot2(z1, main = colnames(Z)[i]) } } + +#' +outlier.isoforest_zscore <- function(X, scale=TRUE, ndim=2, ntrees=10000) { + cX <- X - rowMeans(X, na.rm=TRUE) + cX <- cX[complete.cases(cX),,drop=FALSE] + res <- irlba::irlba(cX, nv=3) + V <- res$v + ndim <- min(ndim, min(dim(cX))) + model <- isotree::isolation.forest(V, ndim=ndim, ntrees=ntrees, nthreads=8) + scores <- predict(model, V) + scores <- abs(scores - mean(scores, na.rm=TRUE)) + scores <- scores / sd(scores, na.rm=TRUE) + scores +} From a8fc4cf7fe133a6a9f430f89a71d44002d8d8318 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" Date: Thu, 6 Aug 2026 09:06:11 +0000 Subject: [PATCH 2/7] chore: update DESCRIPTION imports [skip ci] --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index deece11d..e066510e 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -69,13 +69,13 @@ Imports: GSEABase, GSVA, harmony, - HiClimR, htmlwidgets, igraph, iheatmapr, IntNMF, IRanges, irlba, + isotree, isva, jsonlite, karyoploteR, From c4d2fb8027f1af880869bceb21d19156e35b55b2 Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" Date: Thu, 6 Aug 2026 11:42:54 +0000 Subject: [PATCH 3/7] chore: update documentation --- man/getExampleFeatures.Rd | 2 +- man/getOrganismGO.Rd | 3 ++- man/getSpeciesAliases.Rd | 11 +++++++++++ 3 files changed, 14 insertions(+), 2 deletions(-) create mode 100644 man/getSpeciesAliases.Rd diff --git a/man/getExampleFeatures.Rd b/man/getExampleFeatures.Rd index a8342120..98069d32 100644 --- a/man/getExampleFeatures.Rd +++ b/man/getExampleFeatures.Rd @@ -4,7 +4,7 @@ \alias{getExampleFeatures} \title{Return n example features (symbols) for given organism} \usage{ -getExampleFeatures(organism, n, protein.coding = TRUE, type = "SYMBOL") +getExampleFeatures(organism, n = 20, db = c("gprofiler", "orgdb")) } \description{ Return n example features (symbols) for given organism diff --git a/man/getOrganismGO.Rd b/man/getOrganismGO.Rd index 6a4b3e3e..0ed8804d 100644 --- a/man/getOrganismGO.Rd +++ b/man/getOrganismGO.Rd @@ -10,7 +10,8 @@ getOrganismGO( features = NULL, minsize = 3, batch_size = 2000, - db = c("annothub", "gprofiler") + db = c("annothub", "gprofiler"), + include_iea = TRUE ) } \description{ diff --git a/man/getSpeciesAliases.Rd b/man/getSpeciesAliases.Rd new file mode 100644 index 00000000..dfe3b5f0 --- /dev/null +++ b/man/getSpeciesAliases.Rd @@ -0,0 +1,11 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/pgx-annot-utils.R +\name{getSpeciesAliases} +\alias{getSpeciesAliases} +\title{Show all aliases} +\usage{ +getSpeciesAliases(species) +} +\description{ +Show all aliases +} From ae9748f3cc222913628182b7d5c4784c5195a9e6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Xavier=20Escrib=C3=A0=20Montagut?= Date: Mon, 10 Aug 2026 09:54:38 +0200 Subject: [PATCH 4/7] fix outlier detection blockers from PR review - col default was `col=col`, so any plot=TRUE call errored on a missing argument. Default to "grey70", matching plotOutlierScores. - subsetting Z dropped to a vector for single-method calls, breaking rowMeans(). Added drop = FALSE. - gate the isoforest behind "z.isoforest" %in% methods. It ran on every call, including the default where its column is discarded, which made isotree a hard dependency of the default path. Co-Authored-By: Claude Opus 5 --- R/pgx-outlier.R | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/R/pgx-outlier.R b/R/pgx-outlier.R index 95819d3b..47527dcf 100644 --- a/R/pgx-outlier.R +++ b/R/pgx-outlier.R @@ -7,7 +7,7 @@ detectOutlierSamples <- function(X, methods = c("z.correlation", "z.distance", "z.features", "z.isoforest")[1:3], - col=col, plot = TRUE, par = TRUE) { + col = "grey70", plot = TRUE, par = TRUE) { if(is.null(methods)) { methods = c("z.correlation", "z.distance", "z.features", "z.isoforest") @@ -38,12 +38,15 @@ detectOutlierSamples <- function(X, xz <- colMeans(xz, na.rm = TRUE) z3 <- abs(xz - median(xz, na.rm = TRUE)) / mad(xz, na.rm = TRUE) - ## isoforest z-score - z4 <- outlier.isoforest_zscore(X, scale=TRUE, ndim=2, ntrees=10000) - - Z <- cbind(z1, z2, z3, z4) - colnames(Z) <- c("z.correlation", "z.distance", "z.features", "z.isoforest") - Z <- Z[,which(colnames(Z) %in% methods)] + ## isoforest z-score. only on request: fitting 10k trees is expensive + if ("z.isoforest" %in% methods) { + z4 <- outlier.isoforest_zscore(X, scale=TRUE, ndim=2, ntrees=10000) + } + + ## NULL columns are dropped by cbind(), so unrequested methods vanish here + Z <- cbind(z.correlation = z1, z.distance = z2, z.features = z3, + z.isoforest = z4) + Z <- Z[, which(colnames(Z) %in% methods), drop = FALSE] zz <- rowMeans(Z, na.rm = TRUE) z0 <- 0.1 * mean(Z, na.rm = TRUE) From 2252aadfc276b0b5b1fc490eb8d67a7f3d915df7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Xavier=20Escrib=C3=A0=20Montagut?= Date: Mon, 10 Aug 2026 13:08:32 +0200 Subject: [PATCH 5/7] fix review blockers on the new outlier method - expose outlier_methods on pgx.preprocess() so z.isoforest is reachable. Defaults to the current three methods, not NULL: detectOutlierSamples() reads NULL as "all methods", which would enable the isoforest for every existing caller. - validate methods with match.arg(). An unmatched name silently produced a zero-column Z and NaN z-scores (reading as "no outliers") with plot=FALSE, and died in plotOutlierScores() with "subscript out of bounds" with plot=TRUE. - swap irlba for stdlib svd() in outlier.isoforest_zscore. nv=3 was hardcoded, so it aborted at <=3 samples and warned at 4. Callers cap the input at 1000 rows, so a truncated solver gains nothing. - clamp ndim against ncol(V), the matrix it is applied to, not the input dims. - drop the unused scale argument. - seq_len() instead of 1:ncol(Z). - cover all of the above in test-pgx-outlier.R. Co-Authored-By: Claude Opus 5 --- R/pgx-outlier.R | 21 ++++++----- R/pgx-preprocess.R | 8 ++++- tests/testthat/test-pgx-outlier.R | 58 +++++++++++++++++++++++++++++++ 3 files changed, 77 insertions(+), 10 deletions(-) create mode 100644 tests/testthat/test-pgx-outlier.R diff --git a/R/pgx-outlier.R b/R/pgx-outlier.R index 47527dcf..c4fc07aa 100644 --- a/R/pgx-outlier.R +++ b/R/pgx-outlier.R @@ -9,9 +9,11 @@ detectOutlierSamples <- function(X, "z.features", "z.isoforest")[1:3], col = "grey70", plot = TRUE, par = TRUE) { - if(is.null(methods)) { - methods = c("z.correlation", "z.distance", "z.features", "z.isoforest") - } + all.methods <- c("z.correlation", "z.distance", "z.features", "z.isoforest") + if (is.null(methods)) methods <- all.methods + ## errors on a typo. without this an unmatched name silently yields a + ## zero-column Z and NaN z-scores, which read as "no outliers". + methods <- match.arg(methods, all.methods, several.ok = TRUE) ## correlation and distance X <- head(X[order(-matrixStats::rowSds(X, na.rm = TRUE)), ], 1000) @@ -40,7 +42,7 @@ detectOutlierSamples <- function(X, ## isoforest z-score. only on request: fitting 10k trees is expensive if ("z.isoforest" %in% methods) { - z4 <- outlier.isoforest_zscore(X, scale=TRUE, ndim=2, ntrees=10000) + z4 <- outlier.isoforest_zscore(X, ndim=2, ntrees=10000) } ## NULL columns are dropped by cbind(), so unrequested methods vanish here @@ -72,19 +74,20 @@ plotOutlierScores <- function(res.outliers, z.threshold = c(3, 6, 9), } barplot2(zz, main = "z.outlier (mean)") barplot2(zz2, main = "z.outlier (geom.mean)") - for (i in 1:ncol(Z)) { + for (i in seq_len(ncol(Z))) { z1 <- Z[, i] barplot2(z1, main = colnames(Z)[i]) } } #' -outlier.isoforest_zscore <- function(X, scale=TRUE, ndim=2, ntrees=10000) { +outlier.isoforest_zscore <- function(X, ndim=2, ntrees=10000) { cX <- X - rowMeans(X, na.rm=TRUE) cX <- cX[complete.cases(cX),,drop=FALSE] - res <- irlba::irlba(cX, nv=3) - V <- res$v - ndim <- min(ndim, min(dim(cX))) + ## ponytail: plain svd. callers cap cX at 1000 rows so a truncated + ## solver buys nothing, and irlba(nv=3) aborted below 4 samples. + V <- svd(cX, nu = 0, nv = min(3, ncol(cX)))$v + ndim <- min(ndim, ncol(V)) model <- isotree::isolation.forest(V, ndim=ndim, ntrees=ntrees, nthreads=8) scores <- predict(model, V) scores <- abs(scores - mean(scores, na.rm=TRUE)) diff --git a/R/pgx-preprocess.R b/R/pgx-preprocess.R index b5f81da8..877610a6 100644 --- a/R/pgx-preprocess.R +++ b/R/pgx-preprocess.R @@ -33,6 +33,9 @@ #' \item{impute_method}{Imputation method passed to `imputeMissing`. Default "SVD2".} #' \item{remove_outliers}{Drop outlier samples. Default FALSE.} #' \item{outlier_threshold}{z-score cutoff for `detectOutlierSamples`. Default 3.} +#' \item{outlier_methods}{z-score methods for `detectOutlierSamples`: any of +#' "z.correlation", "z.distance", "z.features", "z.isoforest". Default the +#' first three; "z.isoforest" is opt-in as it fits an isolation forest.} #' \item{meth_type}{Methylation array type for `normalizeMethylation`. Default NULL.} #' } #' @@ -59,6 +62,9 @@ pgx.preprocess <- function(counts, impute_method = "SVD2", remove_outliers = FALSE, outlier_threshold = 3, + ## NB: not NULL. detectOutlierSamples() reads NULL as "all methods", + ## which would switch on the isoforest behind the caller's back. + outlier_methods = c("z.correlation", "z.distance", "z.features"), meth_type = NULL ), options @@ -194,7 +200,7 @@ pgx.preprocess <- function(counts, X <- playbase::imputeMissing(X, method = "SVD2") } } - res <- playbase::detectOutlierSamples(X, plot = FALSE) + res <- playbase::detectOutlierSamples(X, methods = opt$outlier_methods, plot = FALSE) is.outlier <- (res$z.outlier > opt$outlier_threshold) if (any(is.outlier) && !all(is.outlier)) { X <- X[, which(!is.outlier), drop = FALSE] diff --git a/tests/testthat/test-pgx-outlier.R b/tests/testthat/test-pgx-outlier.R new file mode 100644 index 00000000..15226ad2 --- /dev/null +++ b/tests/testthat/test-pgx-outlier.R @@ -0,0 +1,58 @@ +## Covers detectOutlierSamples() method selection and the isoforest path, +## which pgx.preprocess() enables via opt$outlier_methods. + +mk_X <- function(n) { + set.seed(1) + X <- matrix(rnorm(500 * n), 500, n, + dimnames = list(paste0("g", 1:500), paste0("s", 1:n))) + X[, n] <- X[, n] + 4 ## last sample is the outlier + X +} + +test_that("detectOutlierSamples returns one Z column per requested method", { + X <- mk_X(8) + res <- playbase::detectOutlierSamples(X, plot = FALSE) + expect_equal(colnames(res$Z), c("z.correlation", "z.distance", "z.features")) + expect_equal(unname(which.max(res$z.outlier)), 8L) + + ## a single method must stay a matrix, else rowMeans() fails + res1 <- playbase::detectOutlierSamples(X, methods = "z.distance", plot = FALSE) + expect_true(is.matrix(res1$Z)) + expect_equal(colnames(res1$Z), "z.distance") + expect_length(res1$z.outlier, ncol(X)) +}) + +test_that("detectOutlierSamples rejects an unknown method", { + ## used to return NaN z-scores silently, reading as "no outliers" + expect_error( + playbase::detectOutlierSamples(mk_X(8), methods = "z.correlaton", plot = FALSE), + "should be one of" + ) +}) + +test_that("isoforest is off by default and scores the outlier when asked", { + skip_if_not_installed("isotree") + X <- mk_X(8) + expect_false("z.isoforest" %in% colnames(playbase::detectOutlierSamples(X, plot = FALSE)$Z)) + + set.seed(42) + res <- playbase::detectOutlierSamples(X, methods = "z.isoforest", plot = FALSE) + expect_equal(colnames(res$Z), "z.isoforest") + expect_equal(unname(which.max(res$z.outlier)), 8L) + + ## used to abort in irlba(nv = 3) with fewer than 4 samples + small <- playbase::detectOutlierSamples(mk_X(3), methods = "z.isoforest", plot = FALSE) + expect_length(small$z.outlier, 3L) +}) + +test_that("pgx.preprocess plumbs outlier_methods to detectOutlierSamples", { + X <- mk_X(8) + counts <- 2^X + samples <- data.frame(group = rep(c("a", "b"), each = 4), row.names = colnames(X)) + run <- function(...) playbase::pgx.preprocess(counts, samples, contrasts = NULL, + options = list(remove_outliers = TRUE, outlier_threshold = 3, ...)) + + ## the option must reach detectOutlierSamples()' validation, not be ignored + expect_error(run(outlier_methods = "nonsense"), "should be one of") + expect_no_error(run(outlier_methods = "z.distance")) +}) From 9c0fbdc8b7e72dd46be9d48793d23f5d558681f7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Xavier=20Escrib=C3=A0=20Montagut?= Date: Mon, 10 Aug 2026 13:38:46 +0200 Subject: [PATCH 6/7] fix NaN z-scores and plotting on degenerate input Found by a second review pass over the outlier changes. - mad() is 0 as soon as over half the values tie, e.g. the same sample uploaded twice. All three z-scores then evaluated 0/0, so z.outlier was NaN and pgx.preprocess died on "missing value where TRUE/FALSE needed". The three methods shared one abs(x-median)/mad expression, so the guard lives in one local zscore(): fall back to the mean absolute deviation, and to zeros only when every value is identical. Returning zeros directly would instead mask a real outlier whenever over half the samples tie. Bit-identical whenever mad() > 0. - make is.outlier NA-safe so no non-finite score can reach that if() again, whatever produces it. - plotOutlierScores: compute ylim from finite values only (non-finite Z aborted with "need finite 'ylim' values", reachable since par gained a working default), and restore the caller's par() on exit. - test the above, plus the tied-sample regression and plot=TRUE on degenerate input. Drop the isoforest which.max assertion: it was seed-luck, and the score's scale and sign are still under review. Co-Authored-By: Claude Opus 5 --- R/pgx-outlier.R | 31 +++++++++++++++---- R/pgx-preprocess.R | 4 ++- tests/testthat/test-pgx-outlier.R | 50 +++++++++++++++++++++++++++++-- 3 files changed, 75 insertions(+), 10 deletions(-) diff --git a/R/pgx-outlier.R b/R/pgx-outlier.R index c4fc07aa..74202088 100644 --- a/R/pgx-outlier.R +++ b/R/pgx-outlier.R @@ -23,22 +23,35 @@ detectOutlierSamples <- function(X, distX <- as.matrix(dist(t(X))) z1=z2=z3=z4=NULL - + + ## robust z-score. mad() is 0 as soon as over half the values tie (e.g. a + ## sample uploaded twice), which used to give 0/0 = NaN for every sample and + ## crash the caller. Fall back to the mean absolute deviation, and to zeros + ## when every value really is identical (nothing deviates, so nothing is an + ## outlier). Unchanged whenever mad() > 0. + zscore <- function(x) { + dx <- abs(x - median(x, na.rm = TRUE)) + s <- mad(x, na.rm = TRUE) + if (s == 0) s <- mean(dx, na.rm = TRUE) / 0.7979 ## consistency constant + if (s == 0) return(dx * 0) + dx / s + } + ## z-score based on correlation cor.median <- apply(abs(corX), 1, median, na.rm = TRUE) x1 <- (cor.median - mean(cor.median, na.rm = TRUE)) - z1 <- abs(x1 - median(x1, na.rm = TRUE)) / mad(x1, na.rm = TRUE) + z1 <- zscore(x1) ## z-score based on euclidean distance dist.max <- apply(distX, 1, max, na.rm = TRUE) dist.q10 <- apply(distX, 1, quantile, probs = 0.1, na.rm = TRUE) dist.r <- dist.q10 / dist.max - z2 <- abs(dist.r - median(dist.r, na.rm = TRUE)) / mad(dist.r, na.rm = TRUE) + z2 <- zscore(dist.r) ## gene-wise z-score xz <- abs(X - rowMeans(X, na.rm = TRUE)) / matrixStats::rowSds(X, na.rm = TRUE) xz <- colMeans(xz, na.rm = TRUE) - z3 <- abs(xz - median(xz, na.rm = TRUE)) / mad(xz, na.rm = TRUE) + z3 <- zscore(xz) ## isoforest z-score. only on request: fitting 10k trees is expensive if ("z.isoforest" %in% methods) { @@ -62,13 +75,19 @@ detectOutlierSamples <- function(X, #' @export plotOutlierScores <- function(res.outliers, z.threshold = c(3, 6, 9), col = "grey70", par = TRUE) { - if (par) par(mfrow = c(2, 3), mar = c(8, 4, 2, 2)) + if (par) { + ## restore the caller's layout, we are a guest in their device + opar <- graphics::par(mfrow = c(2, 3), mar = c(8, 4, 2, 2)) + on.exit(graphics::par(opar)) + } Z <- res.outliers$Z zz <- res.outliers$z.outlier zz2 <- res.outliers$z.outlier2 barplot2 <- function(x, ...) { barplot(x, col = col, las = 3, - ylim = c(0, max(10, max(Z))), + ## finite only: an all-NA/Inf column gave "need finite 'ylim' values". + ## max(10, numeric(0)) is 10, so a fully non-finite Z still plots. + ylim = c(0, max(10, Z[is.finite(Z)])), ylab = "z-score", ...) abline(h = z.threshold, lty = 3, col = "red") } diff --git a/R/pgx-preprocess.R b/R/pgx-preprocess.R index 877610a6..2ce9eb07 100644 --- a/R/pgx-preprocess.R +++ b/R/pgx-preprocess.R @@ -201,7 +201,9 @@ pgx.preprocess <- function(counts, } } res <- playbase::detectOutlierSamples(X, methods = opt$outlier_methods, plot = FALSE) - is.outlier <- (res$z.outlier > opt$outlier_threshold) + ## NA-safe: a non-finite score must not be read as an outlier, and must not + ## reach the if() below as NA ("missing value where TRUE/FALSE needed"). + is.outlier <- !is.na(res$z.outlier) & (res$z.outlier > opt$outlier_threshold) if (any(is.outlier) && !all(is.outlier)) { X <- X[, which(!is.outlier), drop = FALSE] counts <- counts[, colnames(X), drop = FALSE] diff --git a/tests/testthat/test-pgx-outlier.R b/tests/testthat/test-pgx-outlier.R index 15226ad2..18b6cb96 100644 --- a/tests/testthat/test-pgx-outlier.R +++ b/tests/testthat/test-pgx-outlier.R @@ -30,21 +30,65 @@ test_that("detectOutlierSamples rejects an unknown method", { ) }) -test_that("isoforest is off by default and scores the outlier when asked", { +test_that("isoforest is off by default and runs when asked", { skip_if_not_installed("isotree") X <- mk_X(8) expect_false("z.isoforest" %in% colnames(playbase::detectOutlierSamples(X, plot = FALSE)$Z)) - set.seed(42) + ## NB: structure only. The isoforest's score scale and sign are under review + ## (it is two-sided and sd-normalized, unlike the mad-based methods), so this + ## deliberately does not assert which sample scores highest. res <- playbase::detectOutlierSamples(X, methods = "z.isoforest", plot = FALSE) expect_equal(colnames(res$Z), "z.isoforest") - expect_equal(unname(which.max(res$z.outlier)), 8L) + expect_true(all(is.finite(res$z.outlier))) ## used to abort in irlba(nv = 3) with fewer than 4 samples small <- playbase::detectOutlierSamples(mk_X(3), methods = "z.isoforest", plot = FALSE) expect_length(small$z.outlier, 3L) }) +test_that("the mad-based methods score the injected outlier well above threshold", { + res <- playbase::detectOutlierSamples(mk_X(8), plot = FALSE) + expect_gt(res$z.outlier[["s8"]], 3) + expect_true(all(res$z.outlier[paste0("s", 1:7)] < 3)) +}) + +test_that("tied samples give finite scores instead of NaN", { + ## mad() is 0 when over half the values tie, which used to yield 0/0 for + ## every method and crash pgx.preprocess with "missing value where + ## TRUE/FALSE needed". Same sample uploaded 3x, twice. + set.seed(1) + X <- matrix(rnorm(500 * 6), 500, 6, + dimnames = list(paste0("g", 1:500), paste0("s", 1:6))) + X[, 2] <- X[, 1]; X[, 3] <- X[, 1] + X[, 5] <- X[, 4]; X[, 6] <- X[, 4] + + res <- playbase::detectOutlierSamples(X, plot = FALSE) + expect_true(all(is.finite(res$Z))) + expect_true(all(is.finite(res$z.outlier))) + + samples <- data.frame(group = rep(c("a", "b"), each = 3), row.names = colnames(X)) + out <- playbase::pgx.preprocess(2^X, samples, contrasts = NULL, + options = list(remove_outliers = TRUE, outlier_threshold = 3)) + expect_equal(ncol(out$X), 6L) ## identical samples, none is an outlier +}) + +test_that("plotting works on degenerate input and restores par", { + pdf(NULL) + on.exit({ dev.off(); unlink("Rplots.pdf") }, add = TRUE) + before <- graphics::par("mfrow") + + ## non-finite Z used to abort with "need finite 'ylim' values" + res <- playbase::detectOutlierSamples(mk_X(3), plot = TRUE) + expect_length(res$z.outlier, 3L) + expect_equal(graphics::par("mfrow"), before) + + ## all four methods, 6 panels + skip_if_not_installed("isotree") + expect_length(playbase::detectOutlierSamples(mk_X(8), methods = NULL, plot = TRUE)$z.outlier, 8L) + expect_equal(graphics::par("mfrow"), before) +}) + test_that("pgx.preprocess plumbs outlier_methods to detectOutlierSamples", { X <- mk_X(8) counts <- 2^X From 7c6c90f6b7288e07b3dfe2cde15b489713fa60ef Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" Date: Mon, 10 Aug 2026 12:11:26 +0000 Subject: [PATCH 7/7] chore: update documentation --- man/pgx.preprocess.Rd | 3 +++ 1 file changed, 3 insertions(+) diff --git a/man/pgx.preprocess.Rd b/man/pgx.preprocess.Rd index 95903089..0505bf52 100644 --- a/man/pgx.preprocess.Rd +++ b/man/pgx.preprocess.Rd @@ -35,6 +35,9 @@ pgx.preprocess( \item{impute_method}{Imputation method passed to `imputeMissing`. Default "SVD2".} \item{remove_outliers}{Drop outlier samples. Default FALSE.} \item{outlier_threshold}{z-score cutoff for `detectOutlierSamples`. Default 3.} + \item{outlier_methods}{z-score methods for `detectOutlierSamples`: any of + "z.correlation", "z.distance", "z.features", "z.isoforest". Default the + first three; "z.isoforest" is opt-in as it fits an isolation forest.} \item{meth_type}{Methylation array type for `normalizeMethylation`. Default NULL.} }} }