From 9f9ae121d831f474b8a270dcfa8e532ca1ee96c9 Mon Sep 17 00:00:00 2001 From: aviezerl Date: Wed, 15 Jul 2026 15:48:45 +0300 Subject: [PATCH] fix: error on unequal-length sequences in calc_seq_pwm calc_seq_pwm() (and extract_pwm(), which wraps it) built a single rectangular one-hot matrix via do.call(rbind, ...). With sequences of different lengths rbind recycled the shorter rows up to the longest length, so the shorter sequences got wrong, batch-dependent scores - surfaced only as an easy-to-miss rbind warning. Fail loudly with a clear message pointing to compute_pwm() (which handles variable lengths correctly) instead of returning garbage. Claude-Session: https://claude.ai/code/session_01PK3qefBGDoBwd9w5226FEq --- NEWS.md | 5 +++++ R/pwm.R | 7 +++++++ tests/testthat/test-pwm.R | 16 ++++++++++++++++ 3 files changed, 28 insertions(+) diff --git a/NEWS.md b/NEWS.md index b66c487..6c8813c 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,10 @@ # prego 0.0.10 +* Fix: `calc_seq_pwm()` / `extract_pwm()` now error clearly when given sequences + of unequal length instead of silently recycling the shorter ones (via `rbind`) + and returning wrong scores. These functions build a single rectangular one-hot + matrix and require equal-length sequences; for variable lengths use + `compute_pwm()`. * Fix: PWM scoring no longer opens thousands of threads / fails on core-limited machines. `compute_pwm()`, `compute_local_pwm()` and `calc_seq_pwm()` (hence `extract_pwm()`) run a small per-sequence BLAS `dgemm` inside an diff --git a/R/pwm.R b/R/pwm.R index d1ad9a1..3655b33 100644 --- a/R/pwm.R +++ b/R/pwm.R @@ -48,6 +48,13 @@ calc_seq_pwm <- function(sequences, mdb, bidirect = TRUE) { stop("sequences must be a character vector") } + # calc_seq_pwm builds a single rectangular one-hot matrix, so all sequences + # must be the same length. Unequal lengths would silently recycle the shorter + # rows (via rbind) and return wrong scores - fail loudly instead. + if (length(unique(nchar(sequences))) > 1) { + cli_abort("All {.field sequences} must have the same length for {.fn extract_pwm}. For sequences of different lengths, use {.fn compute_pwm}.") + } + # Convert sequences to uppercase and save original names sequences <- toupper(sequences) seq_names <- names(sequences) diff --git a/tests/testthat/test-pwm.R b/tests/testthat/test-pwm.R index 20d2b97..a6b3214 100644 --- a/tests/testthat/test-pwm.R +++ b/tests/testthat/test-pwm.R @@ -33,6 +33,22 @@ test_that("compute_pwm and calc_seq_pwm produce the same results", { expect_true(abs(scores_new - scores_old) < 1e-6) }) +test_that("calc_seq_pwm errors on sequences of unequal length", { + test_motif_db <- data.frame( + motif = "test_motif", pos = 1:4, + A = c(0.7, 0.1, 0.1, 0.1), C = c(0.1, 0.7, 0.1, 0.1), + G = c(0.1, 0.1, 0.7, 0.1), T = c(0.1, 0.1, 0.1, 0.7) + ) + test_mdb <- create_motif_db(test_motif_db) + # Unequal lengths used to silently recycle (rbind) and return wrong scores. + expect_error( + calc_seq_pwm(c("ACGTACGT", "ACGTAC"), test_mdb), + "same length" + ) + # Equal lengths still work. + expect_silent(calc_seq_pwm(c("ACGTACGT", "TGCATGCA"), test_mdb)) +}) + test_that("compute_pwm and calc_seq_pwm produce the same results with spatial", { res <- regress_pwm(cluster_sequences_example, cluster_mat_example[, 1], final_metric = "ks", spat_bin_size = 40,