kfukuta/wf-transcriptomes-v202
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Annotation quant issues test [CW-7261]

Browse Source
This commit is contained in:
Sarah Griffiths 2026-05-22 09:38:33 +00:00
parent 18e0c5e92f
commit 7b94eb2ca0
3 changed files with 159 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
.gitlab-ci.yml
View File

@ -67,6 +67,7 @@ docker-run:
- MATRIX_NAME: [
"int_discover_dna", "int_fixed_rna", "int_de_control_vs_control",
"smoke_discover", "smoke_fixed", "smoke_direct_rna", "smoke_de",
"mouse_de_0countquant"
]
rules:
# NOTE As we're overriding the rules block for the included docker-run
@ -170,3 +171,11 @@ docker-run:
test -f ${CI_PROJECT_NAME}/de_analysis/condition_treated_vs_control/results_dge.tsv &&
test -f ${CI_PROJECT_NAME}/de_analysis/condition_treated_vs_control/results_dtu_transcript.tsv &&
[ "$(find ${CI_PROJECT_NAME}/samples -type f -name 'gene_counts.tsv' | wc -l)" -eq 4 ]
# Tests a common error when there are 0 annotation counts for a chunk
- if: $MATRIX_NAME == "mouse_de_0countquant"
variables:
NF_BEFORE_SCRIPT: "mkdir -p ${CI_PROJECT_NAME}/data/ && wget https://ont-exd-int-s3-euwst1-epi2me-labs.s3.amazonaws.com/wf-transcriptomes/mouse_subset_test.tar.gz -O ${CI_PROJECT_NAME}/data/mouse_subset_test.tar.gz && tar -xzvf ${CI_PROJECT_NAME}/data/mouse_subset_test.tar.gz -C ${CI_PROJECT_NAME}/data/"
NF_WORKFLOW_OPTS: "--fastq ${CI_PROJECT_NAME}/data/mouse_subset_test/samples \
--de_analysis --ref_genome ${CI_PROJECT_NAME}/data/mouse_subset_test/mouse_subset.fa \
--ref_annotation ${CI_PROJECT_NAME}/data/mouse_subset_test/mouse_subset.gtf.gz \
--direct_rna --sample_sheet ${CI_PROJECT_NAME}/data/mouse_subset_test/sample_sheet.csv"

43
bin/workflow_glue_r/R/bambu.R
View File

@ -342,6 +342,17 @@ bambu_chunk_rc_files <- function(rc_files, aliases, sample_df) {
}, seqnames, chunk_ids)
}
bambu_known_quant_edge_error_kind <- function(msg) {
if (grepl("Error in filter(., (uniqueStartLengthQuery <= primarySecondaryDistStartEnd", msg, fixed = TRUE)) {
return("uniqueStartLengthQuery_filter")
}
if (grepl("eqClassById` with `y$eqClassById` due to incompatible types.", msg, fixed = TRUE)) {
return("eqClassById_incompatible_types")
}
NA_character_
}
bambu_write_discovery_outputs <- function(out_dir, rc_files, discovered_annotations, chunk_bundles, sample_df) {
saveRDS(rc_files, file.path(out_dir, "bambu_rcfiles.rds"))
saveRDS(discovered_annotations, file.path(out_dir, "bambu_discovered_annotations.rds"))
@ -440,6 +451,7 @@ bambu_run_discover_mode <- function(args, analysis_fn, prepare_annotations_fn, b
)
}
bambu_run_quant_mode <- function(args, analysis_fn) {
chunk_bundle <- readRDS(args$chunk_rds)
discovered_annotations <- readRDS(args$discovered_annotation_rds)
@ -470,7 +482,11 @@ bambu_run_quant_mode <- function(args, analysis_fn) {
sample_names = sample_names
)
} else {
se <- bambu_call_analysis(
# We attempted pre-quant validation (see CW-7261) but we do not seem to be able
# to reliably eliminate these chunk-specific empty/typing edge cases from the discover side.
se <- tryCatch(
bambu_call_analysis(
analysis_fn,
bambu_build_args(
args,
@ -479,7 +495,28 @@ bambu_run_quant_mode <- function(args, analysis_fn) {
discovery = FALSE,
quant = TRUE
)
)
),
error = function(e) {
msg <- conditionMessage(e)
error_kind <- bambu_known_quant_edge_error_kind(msg)
if (is.na(error_kind)) {
stop(e)
}
warning(
sprintf(
"Skipping bambu quant for chunk '%s' (%s): known bambu chunk edge case (%s).",
chunk_bundle$chunk_id,
chunk_bundle$seqname,
error_kind
),
call. = FALSE
)
bambu_empty_quant_se(
discovered_annotations = discovered_annotations,
sample_names = sample_names
)
}
)
}
if (length(sample_names) > 0 && ncol(se) == length(sample_names)) {
@ -860,7 +897,7 @@ bambu_add_library_qc_stats <- function(se, qc_stats = list()) {
sprintf(
paste0(
"Large library size variation detected (%.1fx difference).\n",
" Min: %d, Max: %d reads.\n",
" Min: %.0f, Max: %.0f reads.\n",
" CPM normalization may not be appropriate for such variation."
),
lib_size_ratio,

111
bin/workflow_glue_r/tests/testthat/test_bambu.R
View File

@ -563,7 +563,7 @@ testthat::test_that("quant mode skips chunks with no discovered annotations on t
suppressMessages(main_run_bambu(
args,
analysis_fn = fake_analysis
)),
)) ,
"contains no transcripts on this seqname"
)
@ -589,6 +589,115 @@ testthat::test_that("quant mode skips chunks with no discovered annotations on t
testthat::expect_equal(qc$total_genes_before_filter, 0)
})
testthat::test_that("quant mode catches known uniqueStartLengthQuery/min-Inf edge case and writes empty outputs", {
fixture_dir <- tempfile("bambu-quant-known-edge-")
dir.create(fixture_dir)
chunk_bundle <- list(
chunk_id = "chr1",
seqname = "chr1",
aliases = c("sampleA"),
sample_df = data.frame(alias = c("sampleA"), stringsAsFactors = FALSE),
rc_files = list(sampleA = make_test_tx_se(sample_names = "sampleA")),
annotation_tx_count = 1L
)
chunk_rds <- file.path(fixture_dir, "chr1.rds")
saveRDS(chunk_bundle, chunk_rds)
discovered_annotation_rds <- file.path(fixture_dir, "annotations.rds")
saveRDS(make_test_bambu_row_ranges(fixture_dir), discovered_annotation_rds)
analysis_called <- FALSE
fake_analysis <- function(...) {
analysis_called <<- TRUE
stop(
paste(
"Error in filter(., (uniqueStartLengthQuery <= primarySecondaryDistStartEnd & : In argument: `==...`.",
"Caused by warning in `min()`: no non-missing arguments to min; returning Inf"
),
call. = FALSE
)
}
args <- workflow_glue_r_normalise_args(
list(
mode = "quant",
genome = "genome.fa",
out_dir = file.path(fixture_dir, "out"),
chunk_rds = chunk_rds,
discovered_annotation_rds = discovered_annotation_rds,
transcriptome_mode = "discover",
ndr = NULL,
threads = 2
),
bambu_arg_spec()
)
testthat::expect_warning(
suppressMessages(main_run_bambu(args, analysis_fn = fake_analysis)),
"known bambu chunk edge case"
)
testthat::expect_true(analysis_called)
se <- readRDS(file.path(args$out_dir, "bambu_transcripts.rds"))
testthat::expect_equal(nrow(se), 0)
testthat::expect_equal(colnames(se), c("sampleA"))
})
testthat::test_that("quant mode catches known eqClassById incompatible-type edge case and writes empty outputs", {
fixture_dir <- tempfile("bambu-quant-known-eqclass-edge-")
dir.create(fixture_dir)
chunk_bundle <- list(
chunk_id = "chr2",
seqname = "chr2",
aliases = c("sampleA"),
sample_df = data.frame(alias = c("sampleA"), stringsAsFactors = FALSE),
rc_files = list(sampleA = make_test_tx_se(sample_names = "sampleA")),
annotation_tx_count = 1L
)
chunk_rds <- file.path(fixture_dir, "chr2.rds")
saveRDS(chunk_bundle, chunk_rds)
discovered_annotation_rds <- file.path(fixture_dir, "annotations.rds")
saveRDS(make_test_bambu_row_ranges(fixture_dir), discovered_annotation_rds)
analysis_called <- FALSE
fake_analysis <- function(...) {
analysis_called <<- TRUE
stop(
"Can't join `x$eqClassById` with `y$eqClassById` due to incompatible types.",
call. = FALSE
)
}
args <- workflow_glue_r_normalise_args(
list(
mode = "quant",
genome = "genome.fa",
out_dir = file.path(fixture_dir, "out"),
chunk_rds = chunk_rds,
discovered_annotation_rds = discovered_annotation_rds,
transcriptome_mode = "discover",
ndr = NULL,
threads = 2
),
bambu_arg_spec()
)
testthat::expect_warning(
suppressMessages(main_run_bambu(args, analysis_fn = fake_analysis)),
"known bambu chunk edge case"
)
testthat::expect_true(analysis_called)
se <- readRDS(file.path(args$out_dir, "bambu_transcripts.rds"))
testthat::expect_equal(nrow(se), 0)
testthat::expect_equal(colnames(se), c("sampleA"))
})
testthat::test_that("empty mode writes valid empty outputs including bambu rds files", {
fixture_dir <- tempfile("bambu-empty-mode-")
dir.create(fixture_dir)