diff --git a/R/mic.R b/R/mic.R index f431d3392..65c05f435 100755 --- a/R/mic.R +++ b/R/mic.R @@ -56,6 +56,20 @@ COMMON_MIC_VALUES <- c( 1024, 2048, 4096 ) +# Guidelines publish the twofold dilutions below 0.125 mg/L in rounded notations, e.g., EUCAST prints 0.03 for the dilution that this package holds as 0.032, while laboratories may report 0.03, 0.031, 0.03125 or 0.032. These all denote one dilution. This returns the value of that dilution in COMMON_MIC_VALUES, and leaves all other values as they are. +mic_dilution_level <- function(x) { + notation <- c( + 0.0039, 0.0078, 0.015, 0.0156, 0.03, 0.031, 0.0313, 0.03125, 0.06, 0.063, 0.0625, 0.12 + ) + level <- c( + 0.004, 0.008, 0.016, 0.016, 0.032, 0.032, 0.032, 0.032, 0.064, 0.064, 0.064, 0.125 + ) + x <- as.double(x) + hit <- match(x, notation) + x[!is.na(hit)] <- level[hit[!is.na(hit)]] + x +} + #' Transform Input to Minimum Inhibitory Concentrations (MIC) #' #' This transforms vectors to a new class [`mic`], which treats the input as decimal numbers, while maintaining operators (such as ">=") and only allowing valid MIC values known to the field of (medical) microbiology. diff --git a/R/sir.R b/R/sir.R old mode 100755 new mode 100644 index 3dbe7916d..8db1ec1f0 --- a/R/sir.R +++ b/R/sir.R @@ -1906,9 +1906,18 @@ as_sir_method <- function(method_short, ) } - # EUCAST indicates that there is no S category with an off-scale S breakpoint of 0.001 mg/L (stored as 0.0001 in - # clinical_breakpoints), so that wild-type isolates are I; capped MICs such as "<=0.5" must then not be S or NI (#295) + # clinical_breakpoints holds breakpoints as published, so keep these for sir_interpretation_history() (#331) + breakpoint_S_R_published <- paste0(breakpoints_current[, "breakpoint_S", drop = TRUE], "-", breakpoints_current[, "breakpoint_R", drop = TRUE]) + # EUCAST indicates that there is no S category with an off-scale S breakpoint of 0.001 mg/L, so that wild-type isolates are I, also below that value; capped MICs such as "<=0.5" must then not be S or NI (#295) no_S_category <- method == "mic" & isTRUE(breakpoints_current$breakpoint_S <= 0.001) + if (method == "mic") { + # compare dilutions, not their notation, e.g. an MIC of 0.032 and a breakpoint of 0.03 are the same dilution (#331) + breakpoints_current$breakpoint_S <- mic_dilution_level(breakpoints_current$breakpoint_S) + breakpoints_current$breakpoint_R <- mic_dilution_level(breakpoints_current$breakpoint_R) + if (guideline_current %like% "EUCAST" && !custom_breakpoints_set && isTRUE(breakpoints_current$breakpoint_S == 0.001)) { + breakpoints_current$breakpoint_S <- 0.0001 + } + } middle_sir <- if (isTRUE(breakpoints_current$is_SDD)) as.sir("SDD") else as.sir("I") clsi_closed_R <- !custom_breakpoints_set & guideline_current %like% "CLSI" capped_low_not_R <- (as.character(input_original) %like% "^<[0-9]" & as.double(input_clean) <= breakpoints_current$breakpoint_R) | @@ -2097,7 +2106,7 @@ as_sir_method <- function(method_short, guideline = vectorise_log_entry(guideline_current, length(rows)), ref_table = vectorise_log_entry(breakpoints_current[, "ref_tbl", drop = TRUE], length(rows)), uti = vectorise_log_entry(breakpoints_current[, "uti", drop = TRUE], length(rows)), - breakpoint_S_R = vectorise_log_entry(paste0(breakpoints_current[, "breakpoint_S", drop = TRUE], "-", breakpoints_current[, "breakpoint_R", drop = TRUE]), length(rows)), + breakpoint_S_R = vectorise_log_entry(breakpoint_S_R_published, length(rows)), stringsAsFactors = FALSE ) out <- subset(out, !is.na(input_given)) diff --git a/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R b/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R index afb51f65b..fe61c4625 100644 --- a/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R +++ b/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R @@ -469,6 +469,11 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis # ============================================================================== # 4. Cell helpers # ============================================================================== +# The address of a cell in a workbook, e.g. "B12", kept with every parsed breakpoint so that it can be verified against its source +cell_address <- function(r, c) { + paste0(LETTERS[c], r) +} + get_cell <- function(df, r, c) { hit <- df[df$row == r & df$col == c, ] if (nrow(hit) == 0) { @@ -1002,6 +1007,8 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { disk_dose = NA_character_, breakpoint_S = s_val, breakpoint_R = r_val, + source_S = cell_address(r, 2), + source_R = cell_address(r, 3), is_bracketed = is_bracketed_bp(mic_s$base) || is_bracketed_bp(mic_r$base), note = mic_note ) @@ -1027,6 +1034,8 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { disk_dose = disk_dose, breakpoint_S = s_val, breakpoint_R = r_val, + source_S = cell_address(r, 6), + source_R = cell_address(r, 7), is_bracketed = is_bracketed_bp(disk_s$base) || is_bracketed_bp(disk_r$base), note = disk_note ) @@ -1166,6 +1175,7 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { ab = ab_name, disk_dose = NA_character_, breakpoint_S = mic_val, breakpoint_R = mic_val, + source_S = cell_address(od$mic_row, col_idx), source_R = cell_address(od$mic_row, col_idx), is_bracketed = is_bracketed_bp(mic_cell$base), note = mic_note ) @@ -1187,6 +1197,7 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { ab = ab_name, disk_dose = ab_dose, breakpoint_S = disk_val, breakpoint_R = disk_val, + source_S = cell_address(od$disk_row, col_idx), source_R = cell_address(od$disk_row, col_idx), is_bracketed = is_bracketed_bp(disk_cell$base), note = disk_note ) @@ -1344,6 +1355,8 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { disk_dose = NA_character_, breakpoint_S = s_val, breakpoint_R = r_val, + source_S = cell_address(r, sb$s_col), + source_R = cell_address(r, sb$r_col), is_bracketed = is_bracketed_bp(s_cell$base) || is_bracketed_bp(r_cell$base), note = note_text ) @@ -1554,6 +1567,44 @@ if (!is.null(parse_log$disk_dose_inherited)) { print(parse_log$disk_dose_inherited) } +# Every parsed value must be in the workbook cell that it is said to come from. This reads the workbooks a second time with another package than the parser, so that it does not share the parser's assumptions. Note references are superscript text in the same cell, hence the comparison with the start of the cell. +verify_raw_against_workbooks <- function(raw) { + files <- list.files("data-raw", pattern = "^v_[0-9.]+_.*Breakpoint_Tables[.]xlsx$", full.names = TRUE) + file_of_version <- function(version, file_desc) { + hit <- files[grepl(paste0("v_0?", version, "_"), basename(files), fixed = FALSE) & grepl("Antifungal", basename(files), fixed = TRUE) == grepl("Antifungal", file_desc, fixed = TRUE)] + if (length(hit) != 1) stop("no single workbook found for version ", version, " (", file_desc, ")") + hit + } + keys <- distinct(raw, version, file_desc, sheet) + for (i in seq_len(nrow(keys))) { + rows <- raw[raw$version == keys$version[i] & raw$file_desc == keys$file_desc[i] & raw$sheet == keys$sheet[i], , drop = FALSE] + sheet <- suppressMessages(as.data.frame(readxl::read_excel( + file_of_version(keys$version[i], keys$file_desc[i]), + sheet = keys$sheet[i], col_names = FALSE, col_types = "text", .name_repair = "minimal", range = readxl::cell_limits(c(1, 1), c(NA, NA)) + ))) + cell_text <- function(address) { + r <- as.integer(sub("^[A-Z]+", "", address)) + c <- match(sub("[0-9]+$", "", address), LETTERS) + vapply(seq_along(address), function(j) { + if (r[j] > nrow(sheet) || c[j] > ncol(sheet)) NA_character_ else trimws(gsub("\\s+", " ", sheet[r[j], c[j]])) + }, character(1)) + } + for (side in c("S", "R")) { + value <- trimws(rows[[paste0("breakpoint_", side)]]) + has_value <- !is.na(value) & nzchar(value) + text <- cell_text(rows[[paste0("source_", side)]]) + # numbers can be stored as numeric cells, read as text with a different notation + matches <- startsWith(coalesce(text, ""), coalesce(value, "")) | (suppressWarnings(as.numeric(text)) %in% suppressWarnings(as.numeric(value))) + if (any(has_value & !matches)) { + print(cbind(rows[has_value & !matches, c("guideline", "sheet", "mo", "ab", "method", paste0("breakpoint_", side), paste0("source_", side))], in_cell = text[has_value & !matches])) + stop("parsed breakpoint(s) not found in the workbook cell that they come from, see above") + } + } + } + invisible(TRUE) +} +verify_raw_against_workbooks(breakpoints_eucast_raw) + saveRDS(breakpoints_eucast_raw, "data-raw/breakpoints_eucast_raw.rds") write.csv(breakpoints_eucast_raw, "data-raw/breakpoints_eucast_raw.csv", row.names = FALSE, eol = "\n") @@ -2608,32 +2659,18 @@ breakpoints_eucast <- breakpoints_eucast |> # cannot introduce NAs (blocking rows are NA by design) breakpoint_S = as.numeric(breakpoint_S), breakpoint_R = as.numeric(breakpoint_R), - # EUCAST's "arbitrary, off scale" MIC breakpoint S <= 0.001 mg/L (from - # v10.0 on) is meant to categorise wild-type isolates as "I", not "S". - # MICs of 0.0005 mg/L and below do occur in vitro, so, as before, it is - # stored as 0.0001 mg/L to keep such isolates out of "S". - breakpoint_S = if_else(method == "MIC" & !is.na(breakpoint_S) & breakpoint_S == 0.001, 0.0001, breakpoint_S), # same definition as used for the WHONET-based rows uti = !is.na(site) & site %like% "(UTI|urinary|urine)", is_SDD = FALSE # EUCAST has no "susceptible, dose dependent" category (CLSI-only concept) ) |> - # EUCAST prints the twofold dilutions below 0.125 mg/L rounded, e.g. 0.03 and - # 0.06 for 0.03125 and 0.0625, while this package (as.mic(), and the WHONET- - # based breakpoints of all other guidelines) uses 0.032 and 0.064. Since - # as.sir() compares values literally, an MIC of 0.064 would otherwise fall - # above an S breakpoint of 0.06. MIC breakpoints are therefore set to the - # package's dilution level they denote; a value that is not close to any - # dilution level stops the script. The raw data keep EUCAST's notation. + # MIC breakpoints are stored exactly as EUCAST publishes them, e.g. 0.03 and not 0.032, including the off-scale S <= 0.001 mg/L (#331). as.sir() takes care of different notations of the same dilution, see AMR:::mic_dilution_level(). A value that is not close to any dilution level stops the script. mutate(across(c(breakpoint_S, breakpoint_R), function(x) { if (!all(method %in% c("MIC", "DISK"))) stop("unexpected method") is_mic <- method == "MIC" & !is.na(x) - levels <- AMR:::COMMON_MIC_VALUES - nearest <- vapply(x[is_mic], function(z) levels[which.min(abs(log2(z) - log2(levels)))], numeric(1)) - off_level <- abs(log2(x[is_mic]) - log2(nearest)) > 0.15 + off_level <- !AMR:::mic_dilution_level(x[is_mic]) %in% AMR:::COMMON_MIC_VALUES if (any(off_level)) { stop("MIC breakpoint(s) not on a dilution level: ", paste(unique(x[is_mic][off_level]), collapse = ", ")) } - x[is_mic] <- nearest x })) |> # Greek symbols and EM dash symbols are not allowed by CRAN, so replace them with ASCII: diff --git a/tests/testthat/test-data.R b/tests/testthat/test-data.R index c9cc71edd..aa8408914 100644 --- a/tests/testthat/test-data.R +++ b/tests/testthat/test-data.R @@ -582,6 +582,23 @@ test_that("microorganisms.codes: unique codes for existing names", { expect_true(all(is.na(microorganisms.codes$source) | grepl("^WHONET \\(.+\\)$", microorganisms.codes$source))) }) +test_that("clinical_breakpoints: EUCAST MIC breakpoints are stored as published (#331)", { + skip_on_cran() + cb <- as.data.frame(clinical_breakpoints, stringsAsFactors = FALSE) + cb <- cb[cb$type == "human" & cb$method == "MIC" & cb$guideline %like% "EUCAST" & cb$ref_tbl %like% "file: ", , drop = FALSE] + expect_true(nrow(cb) > 0) + # these notations do not occur in the EUCAST Clinical Breakpoint Tables + published_as <- c("0.032" = 0.032, "0.064" = 0.064, "0.0001" = 0.0001) + for (value in names(published_as)) { + wrong <- which(cb$breakpoint_S == published_as[value] | cb$breakpoint_R == published_as[value]) + expect_true(length(wrong) == 0, info = paste("EUCAST MIC breakpoint", value, "in", paste(unique(cb$guideline[wrong]), collapse = ", "))) + } + # but those as published do, e.g., benzylpenicillin for group A streptococci in v16.1 + expect_true(any(cb$breakpoint_S == 0.03, na.rm = TRUE)) + expect_true(any(cb$breakpoint_S == 0.06, na.rm = TRUE)) + expect_true(any(cb$breakpoint_S == 0.001, na.rm = TRUE)) +}) + test_that("clinical_breakpoints: consistent and unambiguous", { skip_on_cran() cb <- clinical_breakpoints diff --git a/tests/testthat/test-sir.R b/tests/testthat/test-sir.R index b06684ad8..571f223e4 100644 --- a/tests/testthat/test-sir.R +++ b/tests/testthat/test-sir.R @@ -638,6 +638,29 @@ test_that("custom reference_data: host = NA acts as host-agnostic fallback", { expect_equal(as.character(result), "R") }) +test_that("MIC breakpoints are compared as dilutions, not as their notation (#331)", { + skip_on_cran() + sir <- function(x, mo, ab, guideline = "EUCAST 2026") { + as.character(suppressMessages(suppressWarnings(as.sir(as.mic(x), mo = mo, ab = ab, guideline = guideline, info = FALSE)))) + } + # the published value is the stored value + pen <- subset(clinical_breakpoints, guideline == "EUCAST 2026" & type == "human" & method == "MIC" & ab == "PEN" & mo == "B_STRPT_GRPA") + expect_identical(unique(pen$breakpoint_S), 0.03) + # all notations of this dilution are S, the next dilution is not + expect_identical(sir(c(0.03, 0.031, 0.03125, 0.032), "Streptococcus pyogenes", "PEN"), rep("S", 4)) + expect_identical(sir(0.064, "Streptococcus pyogenes", "PEN"), "R") + # the same for an R breakpoint, in the notation of the guideline and of the laboratory + expect_identical(mic_dilution_level(c(0.03, 0.032, 0.06, 0.064, 0.12, 0.125, 0.19, 3)), c(0.032, 0.032, 0.064, 0.064, 0.125, 0.125, 0.19, 3)) + # the off-scale S breakpoint is stored as published, and never leads to S + cip <- subset(clinical_breakpoints, guideline == "EUCAST 2026" & type == "human" & method == "MIC" & ab == "CIP" & mo == "B_PSDMN" & is.na(site)) + expect_identical(unique(cip$breakpoint_S), 0.001) + expect_identical(sir(c(0.0005, 0.001, 0.25), "Pseudomonas aeruginosa", "CIP"), c("I", "I", "I")) + # the log shows the published breakpoint + sir_interpretation_history(clean = TRUE) + sir(0.032, "Streptococcus pyogenes", "PEN") + expect_identical(sir_interpretation_history(clean = TRUE)$breakpoint_S_R, "0.03-0.03") +}) + test_that("Morganella is interpreted as reported, not as M. morganii (#330)", { skip_on_cran() cb <- as.data.frame(clinical_breakpoints, stringsAsFactors = FALSE)