(v3.0.1.9153) Speed up mdro(), as.sir(), wisca() and wisca_plot() and the binding of mo vectors, fix PDR in the CMI 2012 guideline of mdro() being assigned to the wrong rows, and limit the stratified WISCA example to the 5 most common species

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: AMR
-Version: 3.0.1.9152
-Date: 2026-10-15
+Version: 3.0.1.9153
+Date: 2026-10-16
 Title: Antimicrobial Resistance Data Analysis
 Description: Functions to simplify and standardise antimicrobial resistance (AMR)
   data analysis and to work with microbial and antimicrobial properties by

NEWS.md

@@ -1,4 +1,4 @@
-# AMR 3.0.1.9152
+# AMR 3.0.1.9153
 
 Planned as v3.1.0, mid October 2026.
 
@@ -80,6 +80,7 @@ Planned as v3.1.0, mid October 2026.
 * Microorganism coercion (`as.mo()`) and properties (`mo_*()`)
   * Considerably faster through one consistent session cache design, which is cleared when lookup tables change, e.g., after `add_custom_microorganisms()` (#313)
   * `mo_*()` functions now retrieve properties of valid `mo` input directly, and only for unique values
+  * Subsetting, combining, and binding rows of `mo` vectors (e.g., in `rbind()`) is now much faster, as codes are checked against a hashed set instead of the whole `microorganisms` data set
   * The mapping of outdated taxonomic names to their current names is now remembered per code
   * Previously coerced input is now recognised again, e.g., for an unknown species of a known genus, instead of being matched again on every call (#313)
   * `mo_uncertainties()` now also lists uncertainties of previously coerced input
@@ -121,6 +122,7 @@ Planned as v3.1.0, mid October 2026.
   * Fixed a bug where both functions failed without the `antimicrobials` argument if `tidyselect` was not available, e.g., for binary installations without `dplyr`
   * Isolates of unknown species are now included in grouped antibiograms and in `wisca()`, as they already were in traditional antibiograms. They were previously excluded silently, with the warning "invalid microorganism code, NA generated". This can change WISCA coverage estimates.
   * Grouped traditional antibiograms now show the (translated) pathogen names, as ungrouped antibiograms do, instead of microorganism codes
+  * `wisca()` and `wisca_plot()` are considerably faster, especially with syndromic groups, by determining intrinsic resistance and pathogen names only once per regimen and pathogen
 * Resistance and susceptibility calculation
   * New argument `guideline` (default: the `AMR_guideline` option, otherwise EUCAST) in `susceptibility()`, `resistance()`, `count_df()`, `proportion_df()`, `sir_df()`, `ggplot_sir()`, `geom_sir()`, `labels_sir_count()`, `antibiogram()`, `wisca()`, and `bug_drug_combinations()` to interpret the 'I' category correctly per guideline. Functions that return results per category also accept `guideline = NULL` to keep S, I, and R separate, while functions that return a single proportion or count give an informative error.
   * When EUCAST is assumed because neither the `guideline` argument nor the `AMR_guideline` option is set, a message explains how to change this, once per session and not in R Markdown or Quarto documents
@@ -145,6 +147,8 @@ Planned as v3.1.0, mid October 2026.
 * MDRO determination (`mdro()`)
   * Resistance can now be inferred from a drug/inhibitor combination when the column of the base drug is absent (e.g., piperacillin inferred from piperacillin/tazobactam). This is controlled with the new argument `infer_from_combinations` (default: `TRUE`) (#209).
   * BRMO classification now includes bacterial complexes (#275)
+  * Fixed a bug in the CMI 2012 guideline where pandrug-resistant (PDR) isolates were classified as extensively drug-resistant (XDR), and another isolate could be marked as PDR instead, if the isolate was not among the first rows of the data. This can change the PDR and XDR results of `guideline = "CMI 2012"` (the default).
+  * Considerably faster, as the classes per isolate are now determined per column instead of per cell
 * `top_n_microorganisms()`
   * New argument `property_for_each` for sub-grouping within the top *n* groups. Rank ordering is enforced, so that only lower taxonomic ranks are allowed.
   * `property = NULL` is now accepted

R/aa_helper_functions.R

@@ -724,6 +724,16 @@ return_after_integrity_check <- function(value, type, check_vector) {
   value
 }
 
+return_after_mo_integrity_check <- function(value) {
+  # same as return_after_integrity_check(), but with a hashed lookup of the codes, see mo_codes_valid()
+  valid <- is.na(value) | mo_codes_valid(value)
+  if (!all(valid)) {
+    warning_("invalid microorganism code, NA generated")
+    value[!valid] <- NA
+  }
+  value
+}
+
 # transforms data set to a tibble with only ASCII values, to comply with CRAN policies
 dataset_UTF8_to_ASCII <- function(df) {
   trans <- function(vect) {

R/aa_helper_pm_functions.R

@@ -354,8 +354,11 @@ pm_group_data_worker <- function(.data, pm_groups) {
   class(res) <- "data.frame"
   nrow_res <- nrow(res)
   rows <- rep(list(NA), nrow_res)
+  # compute the group keys only once, instead of for every group
+  data_keys <- as.character(interaction(.data[, pm_groups]))
+  res_keys <- as.character(interaction(res[, pm_groups]))
   for (i in seq_len(nrow_res)) {
-    rows[[i]] <- which(interaction(.data[, pm_groups]) %in% interaction(res[i, pm_groups]))
+    rows[[i]] <- which(data_keys %in% res_keys[i])
   }
   res$`.rows` <- rows
   res <- res[do.call(order, lapply(pm_groups, function(x) res[, x])), , drop = FALSE]

R/aa_session_cache.R

@@ -101,7 +101,31 @@ reset_ab_cache <- function() {
   AMR_env$ab_previously_coerced <- new_coercion_cache()
 }
 
+mo_codes_valid <- function(x) {
+  # whether codes are in `AMR_env$MO_lookup`; a hashed set avoids scanning the whole lookup table for every check, as `[<-.mo` runs one for every piece that rbind() binds
+  x <- as.character(x)
+  if (length(x) > 1000) {
+    # for long vectors, base matching is faster than looking up every unique value
+    add_MO_lookup_to_AMR_env()
+    return(x %in% as.character(AMR_env$MO_lookup$mo))
+  }
+  if (is.null(AMR_env$MO_codes_hash)) {
+    add_MO_lookup_to_AMR_env()
+    codes <- unique(as.character(AMR_env$MO_lookup$mo))
+    AMR_env$MO_codes_hash <- list2env(stats::setNames(as.list(rep(TRUE, length(codes))), codes), hash = TRUE, size = length(codes))
+  }
+  ux <- unique(x)
+  valid <- vapply(
+    FUN.VALUE = logical(1),
+    ux,
+    function(code) !is.na(code) && nzchar(code) && exists(code, envir = AMR_env$MO_codes_hash, inherits = FALSE),
+    USE.NAMES = FALSE
+  )
+  valid[match(x, ux)]
+}
+
 reset_mo_cache <- function() {
+  AMR_env$MO_codes_hash <- NULL
   AMR_env$MO_index <- NULL
   AMR_env$MO_translated <- NULL
   AMR_env$mo_accepted <- NULL

R/antibiogram.R

@@ -192,8 +192,8 @@
 #'   wisca = TRUE
 #' )
 #'
-#' # stratified by syndrome or clinical group
-#' out <- wisca(example_isolates,
+#' # stratified by syndrome or clinical group, here for the 5 most common species to keep the plot readable
+#' out <- wisca(top_n_microorganisms(example_isolates, n = 5),
 #'   antimicrobials = c("TZP", "TZP+TOB", "TZP+GEN"),
 #'   syndromic_group = "ward"
 #' )
@@ -1230,21 +1230,31 @@ create_wisca_priors <- function(data, sep) {
   beta_prior_beta <- rep(0.5, n_pathogens)
 
   # strongly informative Beta(1, 9999) for intrinsically resistant bug-drug pairs (Bielicki et al., 2016)
-  is_intrinsic <- vapply(
-    FUN.VALUE = logical(1),
-    seq_len(nrow(data)),
-    function(i) {
-      # split by " + ", or wherever `sep` is set to
-      ab_components <- as.ab(trimws(strsplit(as.character(data$ab[i]), trimws(sep), fixed = TRUE)[[1]]))
-      ab_components <- ab_components[!is.na(ab_components)]
-      length(ab_components) > 0 &&
-        all(vapply(
-          FUN.VALUE = logical(1),
-          ab_components,
-          function(ab) any(AMR::intrinsic_resistant$mo == data$mo[i] & AMR::intrinsic_resistant$ab == ab)
-        ))
-    }
-  )
+  # coerce every unique regimen only once and look up the bug-drug pairs as keys, which is much faster than per row
+  if (is.null(AMR_env$intrinsic_resistant_keys)) {
+    # the intrinsic_resistant data set never changes within a session
+    AMR_env$intrinsic_resistant_keys <- paste(as.character(AMR::intrinsic_resistant$mo), as.character(AMR::intrinsic_resistant$ab))
+  }
+  regimens <- unique(as.character(data$ab))
+  regimen_components <- lapply(regimens, function(regimen) {
+    # split by " + ", or wherever `sep` is set to
+    ab_components <- as.ab(trimws(strsplit(regimen, trimws(sep), fixed = TRUE)[[1]]))
+    as.character(ab_components[!is.na(ab_components)])
+  })
+  components <- regimen_components[match(as.character(data$ab), regimens)]
+  # look up all bug-drug pairs at once, and require all components of a regimen to be intrinsically resistant
+  pair_row <- rep(seq_len(nrow(data)), lengths(components))
+  pair_intrinsic <- paste(as.character(data$mo)[pair_row], unlist(components)) %in% AMR_env$intrinsic_resistant_keys
+  is_intrinsic <- lengths(components) > 0 & !seq_len(nrow(data)) %in% pair_row[!pair_intrinsic]
+  for (i in which(is.na(data$mo))) {
+    # keep the original behaviour for a missing microorganism
+    is_intrinsic[i] <- length(components[[i]]) > 0 &&
+      all(vapply(
+        FUN.VALUE = logical(1),
+        components[[i]],
+        function(ab) any(AMR::intrinsic_resistant$mo == data$mo[i] & AMR::intrinsic_resistant$ab == ab)
+      ))
+  }
   beta_prior_alpha[is_intrinsic] <- 1
   beta_prior_beta[is_intrinsic] <- 9999
 
@@ -1635,7 +1645,9 @@ plot_wisca_susceptibility_incidence <- function(wisca_model, sep) {
     )
   })
   df <- do.call(rbind, df_list)
-  df$pathogen <- mo_shortname(df$pathogen, keep_synonyms = TRUE, info = FALSE)
+  # only translate the unique pathogens, as there is a row for every simulation
+  pathogens <- unique(df$pathogen)
+  df$pathogen <- mo_shortname(pathogens, keep_synonyms = TRUE, info = FALSE)[match(df$pathogen, pathogens)]
 
   # order pathogens by median incidence across all regimens
   med_inc <- tapply(df$incidence, df$pathogen, stats::median)

R/mdro.R

@@ -743,6 +743,19 @@ mdro <- function(x = NULL,
   fluoroquinolones <- c(CIP, ENX, FLE, GAT, GEM, GRX, LVX, LOM, MFX, NOR, OFX, PAZ, PEF, PRU, RFL, SPX, TMX, TVA)
   # nolint end
 
+  # helper function returning a logical matrix (rows x cols) of which cells contain any of `values`
+  column_matches <- function(df, rows, cols, values) {
+    matrix(
+      vapply(
+        FUN.VALUE = logical(length(rows)),
+        cols,
+        function(col) as.character(df[rows, col, drop = TRUE]) %in% values
+      ),
+      nrow = length(rows),
+      dimnames = list(NULL, cols)
+    )
+  }
+
   # helper function for editing the table
   trans_tbl <- function(to, rows, cols, any_all, reason = NULL) {
     cols.bak <- cols
@@ -759,19 +772,18 @@ mdro <- function(x = NULL,
         ),
         stringsAsFactors = FALSE
       )
+      # determine the non-susceptible columns per column instead of per cell, which is much faster
+      group_vct <- unique(unname(cols_ab))
+      nonsus <- column_matches(x, rows, group_vct, search_result)
+      current <- x[rows, "all_nonsusceptible_columns", drop = TRUE]
       x[rows, "all_nonsusceptible_columns"] <<- vapply(
         FUN.VALUE = character(1),
-        rows,
-        function(row, group_vct = cols_ab) {
-          cols_nonsus <- vapply(
-            FUN.VALUE = logical(1),
-            x[row, group_vct, drop = FALSE],
-            function(y) y %in% search_result
-          )
+        seq_along(rows),
+        function(i) {
           paste(
             unique(sort(c(
-              unlist(strsplit(x[row, "all_nonsusceptible_columns", drop = TRUE], ", ", fixed = TRUE)),
-              names(cols_nonsus)[cols_nonsus]
+              unlist(strsplit(current[i], ", ", fixed = TRUE)),
+              group_vct[nonsus[i, ]]
             ))),
             collapse = ", "
           )
@@ -827,51 +839,39 @@ mdro <- function(x = NULL,
         ),
         stringsAsFactors = FALSE
       )
+      # determine the matches per column instead of per cell, which is much faster
+      lst_cols <- unique(unname(lst_vector))
+      available <- column_matches(x, rows, lst_cols, VALID_SIR_LEVELS[VALID_SIR_LEVELS != "NI"])
+      nonsus <- column_matches(x, rows, lst_cols, search_result)
+      count_groups <- function(matches) {
+        rowSums(matrix(
+          vapply(
+            FUN.VALUE = logical(length(rows)),
+            lst,
+            function(group) {
+              group <- group[!is.na(group)]
+              rowSums(matches[, lst_cols %in% group, drop = FALSE]) > 0
+            }
+          ),
+          nrow = length(rows)
+        ))
+      }
+
       x[rows, "classes_in_guideline"] <<- length(lst)
-      x[rows, "classes_available"] <<- vapply(
-        FUN.VALUE = double(1),
-        rows,
-        function(row, group_tbl = lst) {
-          sum(vapply(
-            FUN.VALUE = logical(1),
-            group_tbl,
-            function(group) any(unlist(x[row, group[!is.na(group)], drop = TRUE]) %in% VALID_SIR_LEVELS[VALID_SIR_LEVELS != "NI"])
-          ))
-        }
-      )
+      x[rows, "classes_available"] <<- count_groups(available)
 
       if (isTRUE(verbose)) {
         x[rows, "all_nonsusceptible_columns"] <<- vapply(
           FUN.VALUE = character(1),
-          rows,
-          function(row, group_vct = lst_vector) {
-            cols_nonsus <- vapply(FUN.VALUE = logical(1), x[row, group_vct, drop = FALSE], function(y) y %in% search_result)
-            paste(unique(sort(names(cols_nonsus)[cols_nonsus])), collapse = ", ")
-          }
+          seq_along(rows),
+          function(i) paste(sort(lst_cols[nonsus[i, ]]), collapse = ", ")
         )
       }
-      x[rows, "classes_affected"] <<- vapply(
-        FUN.VALUE = double(1),
-        rows,
-        function(row, group_tbl = lst) {
-          sum(
-            vapply(
-              FUN.VALUE = logical(1),
-              group_tbl,
-              function(group) {
-                any(unlist(x[row, group[!is.na(group)], drop = TRUE]) %in% search_result, na.rm = TRUE)
-              }
-            ),
-            na.rm = TRUE
-          )
-        }
-      )
+      x[rows, "classes_affected"] <<- count_groups(nonsus)
       # for PDR; all drugs are R (or I if combine_SI = FALSE)
-      x_transposed <- as.list(as.data.frame(t(x[rows, lst_vector, drop = FALSE]),
-        stringsAsFactors = FALSE
-      ))
-      row_filter <- vapply(FUN.VALUE = logical(1), x_transposed, function(y) all(y %in% search_result, na.rm = TRUE))
-      x[which(row_filter), "classes_affected"] <<- 999
+      row_filter <- rowSums(!nonsus) == 0
+      # `row_filter` is relative to `rows`, so it must be mapped back to the rows of `x`
+      x[rows[row_filter], "classes_affected"] <<- 999
     }
 
     if (isTRUE(info)) {

R/mo.R

@@ -767,7 +767,8 @@ pillar_shaft.mo <- function(x, ...) {
     mo_cols <- NULL
   }
 
-  all_mos <- c(AMR_env$MO_lookup$mo, NA)
+  # as.character() prevents that c.mo() runs its integrity check on the full lookup table
+  all_mos <- c(as.character(AMR_env$MO_lookup$mo), NA)
   if (!all(x %in% all_mos) ||
     (!is.null(df) && !all(unlist(df[, which(mo_cols), drop = FALSE]) %in% all_mos))) {
     # markup old mo codes
@@ -887,7 +888,8 @@ print.mo <- function(x, print.shortnames = FALSE, ...) {
   }
   x <- as.character(x)
   names(x) <- x_names
-  if (!all(x %in% c(AMR_env$MO_lookup$mo, NA))) {
+  # as.character() prevents that c.mo() runs its integrity check on the full lookup table
+  if (!all(as.character(x) %in% c(as.character(AMR_env$MO_lookup$mo), NA))) {
     warning_(
       "Some MO codes are from another AMR package version. ",
       "Please update the MO codes with {.help [{.fun as.mo}](AMR::as.mo)}.",
@@ -921,7 +923,8 @@ summary.mo <- function(object, ...) {
 #' @noRd
 as.data.frame.mo <- function(x, ...) {
   add_MO_lookup_to_AMR_env()
-  if (!all(x %in% c(AMR_env$MO_lookup$mo, NA))) {
+  # as.character() prevents that c.mo() runs its integrity check on the full lookup table
+  if (!all(as.character(x) %in% c(as.character(AMR_env$MO_lookup$mo), NA))) {
     warning_(
       "The data contains old MO codes (from another AMR package version). ",
       "Please update your MO codes with {.help [{.fun as.mo}](AMR::as.mo)}."
@@ -959,7 +962,7 @@ as.data.frame.mo <- function(x, ...) {
   attributes(y) <- attributes(i)
   # must only contain valid MOs
   add_MO_lookup_to_AMR_env()
-  return_after_integrity_check(y, "microorganism code", as.character(AMR_env$MO_lookup$mo))
+  return_after_mo_integrity_check(y)
 }
 #' @method [[<- mo
 #' @export
@@ -969,7 +972,7 @@ as.data.frame.mo <- function(x, ...) {
   attributes(y) <- attributes(i)
   # must only contain valid MOs
   add_MO_lookup_to_AMR_env()
-  return_after_integrity_check(y, "microorganism code", as.character(AMR_env$MO_lookup$mo))
+  return_after_mo_integrity_check(y)
 }
 #' @method c mo
 #' @export
@@ -979,7 +982,7 @@ c.mo <- function(...) {
   y <- NextMethod()
   attributes(y) <- attributes(x)
   add_MO_lookup_to_AMR_env()
-  return_after_integrity_check(y, "microorganism code", as.character(AMR_env$MO_lookup$mo))
+  return_after_mo_integrity_check(y)
 }
 
 #' @method unique mo

R/sir.R

@@ -1178,8 +1178,10 @@ get_guideline <- function(guideline, reference_data) {
     return(guideline)
   }
   guideline_param <- trimws2(toupper(guideline))
-  latest_clsi <- rev(sort(subset(reference_data, guideline %like% "CLSI")$guideline))[1L]
-  latest_eucast <- rev(sort(subset(reference_data, guideline %like% "EUCAST")$guideline))[1L]
+  # only the unique guidelines are needed, which is much faster than subsetting the whole table
+  guidelines_available <- unique(reference_data$guideline)
+  latest_clsi <- rev(sort(guidelines_available[guidelines_available %like% "CLSI"]))[1L]
+  latest_eucast <- rev(sort(guidelines_available[guidelines_available %like% "EUCAST"]))[1L]
   guideline_param[guideline_param == "CLSI"] <- latest_clsi
   guideline_param[guideline_param == "EUCAST"] <- latest_eucast
   # like 'EUCAST2020', should be 'EUCAST 2020'

R/vctrs.R

@@ -160,7 +160,7 @@ vec_cast.character.mo <- function(x, to, ...) {
 }
 vec_cast.mo.character <- function(x, to, ...) {
   add_MO_lookup_to_AMR_env()
-  return_after_integrity_check(x, "microorganism code", as.character(AMR_env$MO_lookup$mo))
+  return_after_mo_integrity_check(x)
 }
 
 # S3: disk ----

man/antibiogram.Rd

@@ -363,8 +363,8 @@ antibiogram(example_isolates,
   wisca = TRUE
 )
 
-# stratified by syndrome or clinical group
-out <- wisca(example_isolates,
+# stratified by syndrome or clinical group, here for the 5 most common species to keep the plot readable
+out <- wisca(top_n_microorganisms(example_isolates, n = 5),
   antimicrobials = c("TZP", "TZP+TOB", "TZP+GEN"),
   syndromic_group = "ward"
 )

tests/testthat/test-mdro.R

@@ -348,3 +348,19 @@ test_that("test-mdro.R", {
     "AMC, AMX (inferred from AMC/AXS), AXS, DAP, LNZ, TEC, VAN"
   )
 })
+
+test_that("test-mdro.R - PDR does not depend on the row position", {
+  # a pandrug-resistant S. aureus must be PDR regardless of the rows before it, which other microorganisms take
+  abx <- c("GEN", "RIF", "CPT", "OXA", "FOX", "CIP", "MFX", "SXT", "FUS", "VAN", "TEC", "TLV", "TGC", "CLI", "DAP", "ERY", "LNZ", "CHL", "FOS", "QDA", "TCY", "DOX", "MNO")
+  pdr_first <- data.frame(mo = as.mo(c("S. aureus", "E. coli", "E. coli", "E. coli")))
+  pdr_last <- data.frame(mo = as.mo(c("E. coli", "E. coli", "E. coli", "S. aureus")))
+  for (ab in abx) {
+    pdr_first[[ab]] <- as.sir(c("R", "S", "S", "S"))
+    pdr_last[[ab]] <- as.sir(c("S", "S", "S", "R"))
+  }
+  out_first <- suppressMessages(suppressWarnings(mdro(pdr_first, guideline = "CMI 2012", info = FALSE)))
+  out_last <- suppressMessages(suppressWarnings(mdro(pdr_last, guideline = "CMI 2012", info = FALSE)))
+  expect_identical(as.character(out_first[1]), "Pandrug-resistant (PDR)")
+  expect_identical(as.character(out_last[4]), "Pandrug-resistant (PDR)")
+  expect_identical(as.character(out_first[2:4]), as.character(out_last[1:3]))
+})