From 0fa321e3efc8e4148d34ba377805cf9307127022 Mon Sep 17 00:00:00 2001 From: Arvind Iyer Date: Wed, 9 Sep 2026 19:25:33 -0400 Subject: [PATCH] Fix crash from factor sample.class and sparse tmb tables (#11) get.blocks() lost sample names when sample.class/alteration.class was a factor, because which() drops names on a named-factor comparison but not a named-character one. The resulting unnamed blocks propagated into 0-column matrices in template.obj.gen()/generateW_block(), which crashed retrieveOutliers() with "replacement has length zero" several steps downstream. Separately, new.AL.general() summed am$tmb tables by raw vector position, so a tmb table covering only a subset of samples (e.g. only samples with >=1 mutation of a given type) got silently recycled instead of matched by sample identity, corrupting TMB totals without erroring. Both are now handled at ingestion in new.AL.general(): factor covariates are auto-coerced to character, and every tmb table is aligned to the full sample set by name with missing samples zero-filled. Unknown sample ids or duplicate entries in a tmb table now raise a clear error instead of being silently mishandled. Co-Authored-By: Claude Sonnet 5 --- DESCRIPTION | 2 +- NEWS.md | 7 ++ R/gam_utils.r | 2 +- R/selectX_create.r | 72 ++++++++++---- R/selectX_plot.R | 2 +- R/selectX_run.R | 2 +- R/selectX_stats.r | 2 +- man/new.AL.general.Rd | 18 ++-- .../testthat/test-new.AL.general-robustness.R | 96 +++++++++++++++++++ 9 files changed, 174 insertions(+), 29 deletions(-) create mode 100644 tests/testthat/test-new.AL.general-robustness.R diff --git a/DESCRIPTION b/DESCRIPTION index 1a11d55..c299730 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: SelectSim Title: Selected Events Linked by Evolutionary Conditions in Cancer -Version: 0.1.6 +Version: 0.1.7 Authors@R: c( person("Arvind", "Iyer", , "ayalurarvind@gmail.com", role = c("aut", "cre", "cph"), comment = c(ORCID = "0000-0002-8247-700X")), diff --git a/NEWS.md b/NEWS.md index 7bbf2cf..d937f31 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,10 @@ +# SelectSim 0.1.7 + +* Fixed `replacement has length zero` crash in `retrieveOutliers()` caused by a factor `sample.class`/`alteration.class` silently losing sample names in `get.blocks()` (`which()` drops names on a named factor but not a named character vector); factor input is now auto-coerced to character (#11). +* Fixed silent mis-alignment of per-sample TMB totals in `new.AL.general()` when `am$tmb` tables only listed a subset of samples (e.g. only samples with a nonzero mutation count for that type); TMB tables are now matched and zero-filled by sample identity instead of summed by raw vector position (#11). +* `new.AL.general()` now errors clearly if an `am$tmb` table references a sample id absent from the corresponding `am$M` matrix, or contains duplicate sample entries. +* Added regression tests covering factor `sample.class`, partial/sparse `tmb` tables, and unknown sample ids in `tmb`. + # SelectSim 0.1.6 * Removed `reshape2` dependency; replaced with base R equivalents (`tapply`, `as.data.frame.table`). diff --git a/R/gam_utils.r b/R/gam_utils.r index caa3b17..6278190 100644 --- a/R/gam_utils.r +++ b/R/gam_utils.r @@ -2,7 +2,7 @@ # Author : Marco Mina , Arvind Iyer # Project : SelectSim # Desc : Functions to process the maf to gam -# Version : 0.1.6 +# Version : 0.1.7 ### diff --git a/R/selectX_create.r b/R/selectX_create.r index b169a0f..5ff662e 100644 --- a/R/selectX_create.r +++ b/R/selectX_create.r @@ -2,7 +2,7 @@ # Author : Arvind Iyer, Miljan Petrovic # Project : SelectSim # Desc : Implementation of SelectSim algorithm -# Version : 0.1.6 +# Version : 0.1.7 ### #' Create an Alteration Landscape (AL) object @@ -14,13 +14,17 @@ #' alteration matrices, each genes x samples) and \code{tmb} (a named list of data #' frames, one per matrix in \code{M}, each with columns \code{sample} and #' \code{mutation}). The names of \code{M} and \code{tmb} must match. All matrices -#' in \code{M} must have identical row and column names. -#' @param feat.covariates Named character vector of alteration-type annotations, one -#' entry per feature (gene). Names must match rownames of the matrices in \code{M}. -#' If \code{NULL}, all features are labelled \code{"MUT"}. -#' @param sample.covariates Named character vector of sample-type annotations, one -#' entry per sample. Names must match colnames of the matrices in \code{M}. If -#' \code{NULL}, all samples are labelled \code{"sample"}. +#' in \code{M} must have identical row and column names. Each \code{tmb} table may +#' be partial (list only the samples with a nonzero \code{mutation} count for that +#' type) — any sample present in \code{M} but absent from a \code{tmb} table is +#' treated as having \code{mutation = 0} for that type. +#' @param feat.covariates Named character (or factor, auto-coerced) vector of +#' alteration-type annotations, one entry per feature (gene). Names must match +#' rownames of the matrices in \code{M}. If \code{NULL}, all features are labelled +#' \code{"MUT"}. +#' @param sample.covariates Named character (or factor, auto-coerced) vector of +#' sample-type annotations, one entry per sample. Names must match colnames of the +#' matrices in \code{M}. If \code{NULL}, all samples are labelled \code{"sample"}. #' @param min.freq Minimum number of samples a gene must be mutated in (strictly #' greater than) to be retained. Features with \code{rowSums <= min.freq} are #' dropped. @@ -96,24 +100,53 @@ new.AL.general <- function(am, al$alterations$alteration.class <- rep("MUT", nrow(al$am[[1]])) names(al$alterations$alteration.class) <- rownames(al$am[[1]]) } else { + if (is.factor(feat.covariates)) { + if (verbose) message("feat.covariates was a factor; converting to character.") + feat.covariates <- setNames(as.character(feat.covariates), names(feat.covariates)) + } al$alterations$alteration.class <- feat.covariates } if (is.null(sample.covariates)) { al$samples$sample.class <- rep("sample", ncol(al$am[[1]])) names(al$samples$sample.class) <- colnames(al$am[[1]]) } else { + if (is.factor(sample.covariates)) { + if (verbose) message("sample.covariates was a factor; converting to character.") + sample.covariates <- setNames(as.character(sample.covariates), names(sample.covariates)) + } al$samples$sample.class <- sample.covariates } - # set the tumor mutation burden vector + # set the tumor mutation burden vector, aligning every tmb table to col.order + # by sample identity (rather than assuming row order/length already match) and + # zero-filling any sample missing from a given tmb table. al$tmb <- list() for (i in names(am$M)) { - al$tmb[[i]] <- am$tmb[[i]] - } - al$tmb[["total"]] <- c(rep(0, ncol(am$M[[1]]))) - for (i in names(am$tmb)) { - al$tmb[["total"]] <- al$tmb[["total"]] + am$tmb[[i]][, c("mutation")] + tb <- am$tmb[[i]] + tb_sample <- as.character(tb$sample) + unknown <- setdiff(tb_sample, col.order) + if (length(unknown) > 0) { + stop(sprintf( + "am$tmb[['%s']] contains sample id(s) not present in am$M[['%s']]: %s", + i, i, paste(unknown[seq_len(min(5, length(unknown)))], collapse = ", ") + )) + } + if (anyDuplicated(tb_sample) > 0) { + stop(sprintf("am$tmb[['%s']] has duplicate 'sample' entries.", i)) + } + missing <- setdiff(col.order, tb_sample) + if (verbose && length(missing) > 0) { + message(sprintf( + "am$tmb[['%s']]: %d of %d samples had no entry - filled with mutation = 0", + i, length(missing), length(col.order) + )) + } + mutation <- setNames(rep(0, length(col.order)), col.order) + mutation[tb_sample] <- tb$mutation + al$tmb[[i]] <- data.frame(sample = col.order, mutation = mutation, + row.names = col.order, stringsAsFactors = FALSE) } - names(al$tmb$total) <- am$tmb[[1]]$sample + al$tmb[["total"]] <- Reduce(`+`, lapply(al$tmb, function(t) t[col.order, "mutation"])) + names(al$tmb$total) <- col.order class(al) <- "AL" return(al) @@ -144,8 +177,13 @@ get.blocks <- function(al) { al$samples$sample.class <- rep("sample", ncol(al$am[[1]])) names(al$samples$sample.class) <- colnames(al$am[[1]]) } - alteration.class <- al$alterations$alteration.class[rownames(al$am$full)] - sample.class <- al$samples$sample.class[colnames(al$am$full)] + # Coerce defensively to character: which() silently drops names when comparing + # a named factor (but not a named character vector), which would otherwise + # produce unnamed sample/feature blocks downstream. + alteration.class <- as.character(al$alterations$alteration.class[rownames(al$am$full)]) + names(alteration.class) <- rownames(al$am$full) + sample.class <- as.character(al$samples$sample.class[colnames(al$am$full)]) + names(sample.class) <- colnames(al$am$full) feature.blocks <- lapply(unique(alteration.class), function(x) which(alteration.class == x)) names(feature.blocks) <- unique(alteration.class) sample.blocks <- lapply(unique(sample.class), function(x) which(sample.class == x)) diff --git a/R/selectX_plot.R b/R/selectX_plot.R index 4d8b09f..60adde0 100644 --- a/R/selectX_plot.R +++ b/R/selectX_plot.R @@ -2,7 +2,7 @@ # Author : Arvind Iyer, Miljan Petrovic # Project : SelectSim # Desc : The file contains plot related functions -# Version : 0.1.6 +# Version : 0.1.7 ### # Suppress R CMD check notes for ggplot2/ggridges column name variables diff --git a/R/selectX_run.R b/R/selectX_run.R index 14d80d3..1f72214 100644 --- a/R/selectX_run.R +++ b/R/selectX_run.R @@ -2,7 +2,7 @@ # Author : Arvind Iyer, Miljan Petrovic # Project : SelectSim # Desc : Main file which to run the SelecSim algoritm via calling selectX function to create alteration object with background model and funtion to generate the table. -# Version : 0.1.6 +# Version : 0.1.7 # Notes: # - Better Error message and running text # - Edge case: When sample size in less than 2 there is error in computation (need to fix a number to do this analysis) diff --git a/R/selectX_stats.r b/R/selectX_stats.r index 9eb9e59..6f540cc 100644 --- a/R/selectX_stats.r +++ b/R/selectX_stats.r @@ -2,7 +2,7 @@ # Author : Arvind Iyer # Project : SelectSim # Desc : The file which contains the function to generate the stats and table -# Version : 0.1.6 +# Version : 0.1.7 ### #' Initialize an Alteration Landscape Stats (ALS) container diff --git a/man/new.AL.general.Rd b/man/new.AL.general.Rd index aa349ad..443cc84 100644 --- a/man/new.AL.general.Rd +++ b/man/new.AL.general.Rd @@ -17,15 +17,19 @@ new.AL.general( alteration matrices, each genes x samples) and \code{tmb} (a named list of data frames, one per matrix in \code{M}, each with columns \code{sample} and \code{mutation}). The names of \code{M} and \code{tmb} must match. All matrices -in \code{M} must have identical row and column names.} +in \code{M} must have identical row and column names. Each \code{tmb} table may +be partial (list only the samples with a nonzero \code{mutation} count for that +type) — any sample present in \code{M} but absent from a \code{tmb} table is +treated as having \code{mutation = 0} for that type.} -\item{feat.covariates}{Named character vector of alteration-type annotations, one -entry per feature (gene). Names must match rownames of the matrices in \code{M}. -If \code{NULL}, all features are labelled \code{"MUT"}.} +\item{feat.covariates}{Named character (or factor, auto-coerced) vector of +alteration-type annotations, one entry per feature (gene). Names must match +rownames of the matrices in \code{M}. If \code{NULL}, all features are labelled +\code{"MUT"}.} -\item{sample.covariates}{Named character vector of sample-type annotations, one -entry per sample. Names must match colnames of the matrices in \code{M}. If -\code{NULL}, all samples are labelled \code{"sample"}.} +\item{sample.covariates}{Named character (or factor, auto-coerced) vector of +sample-type annotations, one entry per sample. Names must match colnames of the +matrices in \code{M}. If \code{NULL}, all samples are labelled \code{"sample"}.} \item{min.freq}{Minimum number of samples a gene must be mutated in (strictly greater than) to be retained. Features with \code{rowSums <= min.freq} are diff --git a/tests/testthat/test-new.AL.general-robustness.R b/tests/testthat/test-new.AL.general-robustness.R new file mode 100644 index 0000000..e0a73d8 --- /dev/null +++ b/tests/testthat/test-new.AL.general-robustness.R @@ -0,0 +1,96 @@ +make_synthetic_am <- function(n_genes = 20, n_samples = 30, seed = 1) { + set.seed(seed) + samples <- paste0("S", seq_len(n_samples)) + genes <- paste0("G", seq_len(n_genes)) + mk_gam <- function(p) { + matrix(rbinom(n_genes * n_samples, 1, p), nrow = n_genes, + dimnames = list(genes, samples)) + } + M <- list(Nonsense = mk_gam(0.2), Missense = mk_gam(0.1)) + mk_full_tmb <- function(gam) { + csum <- colSums(gam) + data.frame(sample = names(csum), mutation = csum, row.names = names(csum)) + } + list(M = M, tmb = list(Nonsense = mk_full_tmb(M$Nonsense), Missense = mk_full_tmb(M$Missense)), + samples = samples, genes = genes) +} + +test_that("factor sample.class/alteration.class no longer breaks get.blocks() naming", { + fx <- make_synthetic_am() + sample.class.chr <- setNames(rep(c("A", "B"), length.out = length(fx$samples)), fx$samples) + sample.class.factor <- factor(sample.class.chr) + alteration.class <- setNames(rep("MUT", length(fx$genes)), fx$genes) + + al_chr <- new.AL.general(fx, feat.covariates = alteration.class, + sample.covariates = sample.class.chr, min.freq = 1) + al_fac <- new.AL.general(fx, feat.covariates = alteration.class, + sample.covariates = sample.class.factor, min.freq = 1) + + blocks_chr <- get.blocks(al_chr) + blocks_fac <- get.blocks(al_fac) + + # Previously: factor input made which() drop names, leaving blocks unnamed + # (length(names(...)) == 0), which cascaded into empty template matrices. + expect_true(all(vapply(blocks_fac$sample.blocks, function(b) length(names(b)) > 0, logical(1)))) + expect_equal(blocks_chr$sample.blocks, blocks_fac$sample.blocks) + expect_equal(al_chr$tmb$total, al_fac$tmb$total) +}) + +test_that("selectX() runs to completion with a factor sample.class (regression for #11)", { + fx <- make_synthetic_am() + sample.class <- factor(setNames(rep(c("BM", "PCT"), length.out = length(fx$samples)), fx$samples)) + alteration.class <- setNames(rep("MUT", length(fx$genes)), fx$genes) + + expect_no_error({ + result <- selectX(M = fx, sample.class = sample.class, + alteration.class = alteration.class, + min.freq = 1, n.permut = 10, n.cores = 1, verbose = FALSE) + }) + expect_true(is.list(result)) +}) + +test_that("partial (sparse) tmb tables are zero-filled by sample identity, not recycled", { + fx <- make_synthetic_am() + # Drop most Missense rows, keeping only a handful of samples - mirrors the + # reporter's tmb tables, which only listed samples with >=1 mutation of that type. + sparse_tmb <- fx$tmb$Missense[1:5, ] + fx_sparse <- fx + fx_sparse$tmb$Missense <- sparse_tmb + + sample.class <- setNames(rep("sample", length(fx$samples)), fx$samples) + alteration.class <- setNames(rep("MUT", length(fx$genes)), fx$genes) + + expect_no_warning({ + al <- new.AL.general(fx_sparse, feat.covariates = alteration.class, + sample.covariates = sample.class, min.freq = 1) + }) + + # Samples absent from the sparse Missense table should contribute mutation = 0, + # not a positionally-recycled value from an unrelated sample. + dropped_samples <- setdiff(fx$samples, sparse_tmb$sample) + expect_true(all(al$tmb$Missense[dropped_samples, "mutation"] == 0)) + + # Total should equal Nonsense (full) + Missense (zero-filled for dropped + # samples, actual value otherwise), matched by sample identity. + missense_zero_filled <- setNames(rep(0, length(fx$samples)), fx$samples) + missense_zero_filled[sparse_tmb$sample] <- sparse_tmb$mutation + expected_total <- fx$tmb$Nonsense[fx$samples, "mutation"] + missense_zero_filled[fx$samples] + expect_equal(unname(al$tmb$total[fx$samples]), unname(expected_total)) +}) + +test_that("am$tmb with an unknown sample id errors clearly", { + fx <- make_synthetic_am() + bad_tmb <- fx$tmb$Nonsense + bad_tmb$sample[1] <- "NOT_A_REAL_SAMPLE" + rownames(bad_tmb)[1] <- "NOT_A_REAL_SAMPLE" + fx$tmb$Nonsense <- bad_tmb + + sample.class <- setNames(rep("sample", length(fx$samples)), fx$samples) + alteration.class <- setNames(rep("MUT", length(fx$genes)), fx$genes) + + expect_error( + new.AL.general(fx, feat.covariates = alteration.class, + sample.covariates = sample.class, min.freq = 1), + "not present in am\\$M" + ) +})