diff --git a/DESCRIPTION b/DESCRIPTION index d635516f1..a54a14cbb 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,5 +1,5 @@ Package: AMR -Version: 3.0.1.9120 +Version: 3.0.1.9121 Date: 2026-10-12 Title: Antimicrobial Resistance Data Analysis Description: Functions to simplify and standardise antimicrobial resistance (AMR) diff --git a/NEWS.md b/NEWS.md index 932f03b2b..31c9950e4 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,4 +1,4 @@ -# AMR 3.0.1.9120 +# AMR 3.0.1.9121 Planned as v3.1.0, mid October 2026. diff --git a/R/amr_selectors.R b/R/amr_selectors.R index ba65eb3cd..2c6458025 100755 --- a/R/amr_selectors.R +++ b/R/amr_selectors.R @@ -73,7 +73,6 @@ #' # Though they are primarily intended to use for selections and filters. #' # Examples sections below are split into 'dplyr', 'base R', and 'data.table': #' -#' \donttest{ #' \dontrun{ #' # dplyr ------------------------------------------------------------------- #' @@ -94,8 +93,10 @@ #' # get AMR for all aminoglycosides e.g., per ward: #' example_isolates %>% #' group_by(ward) %>% -#' summarise(across(aminoglycosides(), -#' resistance)) +#' summarise(across( +#' aminoglycosides(), +#' resistance +#' )) #' #' # You can combine selectors with '&' to be more specific: #' example_isolates %>% @@ -105,8 +106,10 @@ #' example_isolates %>% #' filter(mo_genus() %in% c("Escherichia", "Klebsiella")) %>% #' group_by(ward) %>% -#' summarise_at(not_intrinsic_resistant(), -#' resistance) +#' summarise_at( +#' not_intrinsic_resistant(), +#' resistance +#' ) #' #' # get susceptibility for antimicrobials whose name contains "trim": #' example_isolates %>% @@ -159,8 +162,9 @@ #' y <- example_isolates %>% filter(carbapenems() == "R") #' z <- example_isolates %>% filter(if_all(carbapenems(), ~ .x == "R")) #' identical(x, y) && identical(y, z) -#' #' } +#' +#' \donttest{ #' # base R ------------------------------------------------------------------ #' #' # select columns 'IPM' (imipenem) and 'MEM' (meropenem) diff --git a/R/antibiogram.R b/R/antibiogram.R index bfe97c2c8..415ba1ab9 100755 --- a/R/antibiogram.R +++ b/R/antibiogram.R @@ -207,7 +207,8 @@ #' top_n_microorganisms(n = 10) %>% #' group_by( #' age_group = age_groups(age, c(25, 50, 75)), -#' gender) %>% +#' gender +#' ) %>% #' wisca(antimicrobials = c("TZP", "TZP+TOB", "TZP+GEN")) #' } #' diff --git a/R/sysdata.rda b/R/sysdata.rda index 887ffd252..b9745cc6b 100755 Binary files a/R/sysdata.rda and b/R/sysdata.rda differ 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 5988101d2..8388548c0 100644 --- a/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R +++ b/data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R @@ -66,18 +66,20 @@ parse_log$no_breakpoints_determined <- NULL # output. split_rich_text <- function(fmt_list, cell_strike = rep(FALSE, length(fmt_list))) { map2_dfr(fmt_list, cell_strike, function(fmt, cell_struck) { - if (is.null(fmt) || nrow(fmt) == 0) + if (is.null(fmt) || nrow(fmt) == 0) { return(tibble(is_rich = FALSE, base_text = NA_character_, note_super = NA_character_, full_text = NA_character_)) + } # A run without its own strike setting inherits the cell's font is_struck <- coalesce(fmt$strike, isTRUE(cell_struck)) fmt <- fmt[!is_struck, , drop = FALSE] - if (nrow(fmt) == 0) + if (nrow(fmt) == 0) { return(tibble(is_rich = TRUE, base_text = NA_character_, note_super = NA_character_, full_text = NA_character_)) + } is_super <- !is.na(fmt$vertAlign) & fmt$vertAlign == "superscript" tibble( is_rich = TRUE, base_text = paste0(fmt$character[!is_super], collapse = ""), - note_super = paste0(fmt$character[is_super], collapse = ""), + note_super = paste0(fmt$character[is_super], collapse = ""), # Text including superscripts, for note bodies, where a superscript is # part of the text rather than a reference, e.g. "Ca2+" or "Mg2+" full_text = paste0(fmt$character, collapse = "") @@ -103,19 +105,20 @@ parse_cells <- function(cells, formats) { mutate( plain_struck = !is_rich & target_strike, base_text = case_when( - plain_struck ~ NA_character_, - is_rich ~ na_if(trimws(base_text), ""), - !is.na(character) ~ trimws(character), - !is.na(numeric) ~ as.character(numeric), - TRUE ~ NA_character_ + plain_struck ~ NA_character_, + is_rich ~ na_if(trimws(base_text), ""), + !is.na(character) ~ trimws(character), + !is.na(numeric) ~ as.character(numeric), + TRUE ~ NA_character_ ), full_text = case_when( - plain_struck ~ NA_character_, - is_rich ~ na_if(trimws(full_text), ""), - TRUE ~ base_text + plain_struck ~ NA_character_, + is_rich ~ na_if(trimws(full_text), ""), + TRUE ~ base_text ), note_super = if_else(is.na(note_super) | note_super == "", - NA_character_, note_super) + NA_character_, note_super + ) ) |> filter(!is.na(base_text) | !is.na(note_super) | !is.na(full_text)) |> select(row, col, base_text, note_super, full_text) @@ -171,12 +174,12 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { } base_bold <- formats$local$font$bold[col_a$local_format_id] base_strike <- coalesce(formats$local$font$strike[col_a$local_format_id], FALSE) - + out <- vector("list", nrow(col_a)) for (i in seq_len(nrow(col_a))) { fmt <- col_a$character_formatted[[i]] row_i <- col_a$row[i] - + if (is.null(fmt) || nrow(fmt) == 0) { # Plain (non-rich-text) cell: nothing to split, whole text is the agent # name, unless the whole cell is struck through @@ -185,19 +188,19 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { } next } - + # Drop struck-through and superscript runs first, exactly as split_rich_text() does is_struck <- coalesce(fmt$strike, base_strike[i]) fmt <- fmt[!is_struck, , drop = FALSE] is_super <- !is.na(fmt$vertAlign) & fmt$vertAlign == "superscript" fmt <- fmt[!is_super, , drop = FALSE] - + if (nrow(fmt) <= 1) { txt <- if (nrow(fmt) == 1) fmt$character[1] else col_a$character[i] out[[i]] <- tibble(row = row_i, ab_text = txt, mo_override = NA_character_) next } - + bold_resolved <- ifelse(is.na(fmt$bold), base_bold[i], fmt$bold) # An indication qualifier may also follow the organism restriction as a @@ -209,7 +212,7 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { trailing_qualifier <- NA_character_ n <- length(bold_resolved) if (n > 2 && bold_resolved[n] && !bold_resolved[n - 1] && - grepl("^\\s*\\([^()]+\\)\\s*$", fmt$character[n])) { + grepl("^\\s*\\([^()]+\\)\\s*$", fmt$character[n])) { trailing_qualifier <- trimws(fmt$character[n]) fmt <- fmt[-n, , drop = FALSE] bold_resolved <- bold_resolved[-n] @@ -223,17 +226,23 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { break } } - + if (is.na(split_at)) { # Either the whole cell is bold (no restriction present -> fine, keep # as a single agent name) or the bold pattern never settles into an # all-non-bold tail, which is the genuinely unexpected case. if (all(bold_resolved) || !any(bold_resolved)) { - out[[i]] <- tibble(row = row_i, - ab_text = paste(c(paste0(fmt$character, collapse = ""), - trailing_qualifier[!is.na(trailing_qualifier)]), - collapse = " "), - mo_override = NA_character_) + out[[i]] <- tibble( + row = row_i, + ab_text = paste( + c( + paste0(fmt$character, collapse = ""), + trailing_qualifier[!is.na(trailing_qualifier)] + ), + collapse = " " + ), + mo_override = NA_character_ + ) next } stop( @@ -244,7 +253,7 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { "trailing organism restriction is not; inspect this cell manually." ) } - + ab_text <- paste0(fmt$character[seq_len(split_at)], collapse = "") if (!is.na(trailing_qualifier)) { ab_text <- paste(gsub(",\\s*$", "", trimws(ab_text)), trailing_qualifier) @@ -252,7 +261,7 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { mo_override <- paste0(fmt$character[(split_at + 1):n], collapse = "") out[[i]] <- tibble(row = row_i, ab_text = ab_text, mo_override = mo_override) } - + bind_rows(out) |> mutate( ab_text = trimws(gsub("[\r\n]+", " ", ab_text)), @@ -281,23 +290,30 @@ split_agent_organism <- function(cells, formats, sheet_name, valid_rows) { # 2. Notes block parser # ============================================================================== parse_notes_block <- function(txt) { - if (is.na(txt) || txt == "") return(list()) - + if (is.na(txt) || txt == "") { + return(list()) + } + txt <- gsub("\r\n", "\n", txt) txt <- gsub("\r", "\n", txt) txt <- gsub("\u00A0", " ", txt) - + find_keys <- function(pattern, txt) { m <- gregexpr(pattern, txt, perl = TRUE)[[1]] - if (m[1] == -1) - return(data.frame(pos = integer(0), key = character(0), - mlen = integer(0), stringsAsFactors = FALSE)) + if (m[1] == -1) { + return(data.frame( + pos = integer(0), key = character(0), + mlen = integer(0), stringsAsFactors = FALSE + )) + } lens <- attr(m, "match.length") keys <- gsub("[.\\s]+$", "", trimws(substring(txt, m, m + lens - 1)), perl = TRUE) - data.frame(pos = as.integer(m), key = keys, mlen = as.integer(lens), - stringsAsFactors = FALSE) + data.frame( + pos = as.integer(m), key = keys, mlen = as.integer(lens), + stringsAsFactors = FALSE + ) } - + # A note key only starts at the beginning of a line or directly after the # end of a sentence, since EUCAST does not always put each note on its own # line (e.g. "...meropenem only. 3.The addition of..." or "...urinary @@ -320,9 +336,11 @@ parse_notes_block <- function(txt) { # Numbered without dot at the start of a line (e.g. C. difficile sheet): # digit, space, uppercase letter r4 <- find_keys("(?:^|(?<=\\n))[ \\t]*\\d+[ \\t]+(?=[A-Z])", txt) - + all_keys <- rbind(r1, r2, r3, r4) - if (nrow(all_keys) == 0) return(list()) + if (nrow(all_keys) == 0) { + return(list()) + } # At the same position, keep the longest (most specific) match, e.g. "1/A." # over "1.", and drop any match that starts within an earlier kept match all_keys <- all_keys[order(all_keys$pos, -all_keys$mlen), ] @@ -333,7 +351,7 @@ parse_notes_block <- function(txt) { keep[seq_len(nrow(all_keys)) > i & all_keys$pos < all_keys$pos[i] + all_keys$mlen[i]] <- FALSE } all_keys <- all_keys[keep, ] - + notes <- list() for (i in seq_len(nrow(all_keys))) { body_start <- all_keys$pos[i] + all_keys$mlen[i] @@ -344,12 +362,16 @@ parse_notes_block <- function(txt) { # reference "A" individually for (k in unique(c(key, if (grepl("/", key)) strsplit(key, "/")[[1]]))) { if (!nzchar(body)) { - stop("parse_notes_block(): note key '", key, "' has an empty body, which means a ", - "cross-reference or other text was taken for a key. Text: '", substr(txt, 1, 300), "'") + stop( + "parse_notes_block(): note key '", key, "' has an empty body, which means a ", + "cross-reference or other text was taken for a key. Text: '", substr(txt, 1, 300), "'" + ) } if (!is.null(notes[[k]]) && !identical(notes[[k]], body)) { - stop("parse_notes_block(): note key '", k, "' occurs more than once with different ", - "bodies. Text: '", substr(txt, 1, 300), "'") + stop( + "parse_notes_block(): note key '", k, "' occurs more than once with different ", + "bodies. Text: '", substr(txt, 1, 300), "'" + ) } notes[[k]] <- body } @@ -366,13 +388,17 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis # Staphylococcus) or "1" with superscript "4.5" (v16.1, H. influenzae) mic_refs <- character(0) disk_refs <- character(0) - if (!is.na(mic_super) && mic_super != "") + if (!is.na(mic_super) && mic_super != "") { mic_refs <- trimws(unlist(strsplit(mic_super, "[,.[:space:]]+"))) - if (!is.na(disk_super) && disk_super != "") + } + if (!is.na(disk_super) && disk_super != "") { disk_refs <- trimws(unlist(strsplit(disk_super, "[,.[:space:]]+"))) + } all_refs <- unique(c(mic_refs, disk_refs)) all_refs <- all_refs[nzchar(all_refs)] - if (length(all_refs) == 0) return(NA_character_) + if (length(all_refs) == 0) { + return(NA_character_) + } # A reference that is not a key of the table's own notes block can be an # abbreviation from the workbook's general "Notes" sheet, e.g. superscript @@ -383,16 +409,18 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis nl <- c(nl, abbreviations[from_abbr]) unresolved <- setdiff(all_refs, names(nl)) if (length(unresolved) > 0) { - stop("resolve_notes(): note reference(s) '", paste(unresolved, collapse = "', '"), - "' not found in the notes block (keys found: ", paste(names(notes_list), collapse = " "), - ") nor in the workbook abbreviations (", paste(names(abbreviations), collapse = " "), ")") + stop( + "resolve_notes(): note reference(s) '", paste(unresolved, collapse = "', '"), + "' not found in the notes block (keys found: ", paste(names(notes_list), collapse = " "), + ") nor in the workbook abbreviations (", paste(names(abbreviations), collapse = " "), ")" + ) } - + # Detect combined keys (e.g. "1/A") combined_keys <- grep("/", names(nl), value = TRUE) used_combined <- character(0) consumed_refs <- character(0) - + for (ck in combined_keys) { parts <- unlist(strsplit(ck, "/")) if (all(parts %in% all_refs)) { @@ -400,12 +428,12 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis consumed_refs <- c(consumed_refs, parts) } } - + remaining_refs <- setdiff(all_refs, consumed_refs) - + parts_out <- character(0) bodies_seen <- character(0) - + # Combined notes for (ck in used_combined) { bodies <- nl[[ck]] @@ -418,7 +446,7 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis } } } - + # Remaining individual notes: prefer duplicated over missing for (ref in remaining_refs) { if (!is.null(nl[[ref]])) { @@ -431,8 +459,10 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis } } } - - if (length(parts_out) == 0) return(NA_character_) + + if (length(parts_out) == 0) { + return(NA_character_) + } paste(parts_out, collapse = " | ") } @@ -441,7 +471,9 @@ resolve_notes <- function(notes_list, mic_super, disk_super, abbreviations = lis # ============================================================================== get_cell <- function(df, r, c) { hit <- df[df$row == r & df$col == c, ] - if (nrow(hit) == 0) return(list(base = NA_character_, super = NA_character_)) + if (nrow(hit) == 0) { + return(list(base = NA_character_, super = NA_character_)) + } base <- hit$base_text[1] super <- hit$note_super[1] # EUCAST occasionally omits the superscript formatting of a note reference @@ -456,7 +488,9 @@ get_cell <- function(df, r, c) { } parse_bp <- function(txt) { - if (is.na(txt) || txt == "") return(NA_character_) + if (is.na(txt) || txt == "") { + return(NA_character_) + } t <- trimws(txt) # placeholders such as "-" (susceptibility testing not recommended), "IE" # or "Note" are kept as they are: the cleanup turns them into blocking rows @@ -485,14 +519,20 @@ is_bracketed_bp <- function(txt) { # marks screening tests in the agent name instead, e.g. "Cefoxitin # (screen only)". The definition is read per workbook by # parse_workbook_notes() and added to the note of every bracketed row. - if (is.na(txt) || txt == "") return(FALSE) + if (is.na(txt) || txt == "") { + return(FALSE) + } grepl("^\\(.+\\)$", trimws(txt)) } format_disk_dose <- function(dose) { - if (is.na(dose) || dose == "") return(NA_character_) + if (is.na(dose) || dose == "") { + return(NA_character_) + } d <- trimws(dose) - if (grepl("mcg|unit", d, ignore.case = TRUE)) return(d) + if (grepl("mcg|unit", d, ignore.case = TRUE)) { + return(d) + } d <- gsub("-", "/", d) paste0(d, " mcg") } @@ -514,12 +554,14 @@ detect_version <- function(cells) { filter(row == 1, !is_blank, col > 1) |> mutate(text = coalesce(character, "")) version_text <- paste(version_cells$text, collapse = " ") - - version <- regmatches(version_text, - regexpr("v\\.?\\s*[\\d.]+", version_text, perl = TRUE)) + + version <- regmatches( + version_text, + regexpr("v\\.?\\s*[\\d.]+", version_text, perl = TRUE) + ) version <- gsub("v\\.?\\s*", "", version) if (length(version) == 0) version <- NA_character_ - + # The guideline year is read directly from the "valid from YYYY-MM-DD" # date in the version string, not derived from the version number: the # two are not reliably related. EUCAST's bacterial clinical breakpoint @@ -539,7 +581,7 @@ detect_version <- function(cells) { year <- regmatches(version_text, regexpr("\\d{4}", version_text)) } guideline <- if (length(year) == 1) paste0("EUCAST ", year) else NA_character_ - + # A human-readable description of the source file/table, for ref_tbl, # e.g. "Clinical Breakpoint Tables v. 16.1" or "Antifungal Clinical # Breakpoint Table v. 12.1" -- everything in the version string between @@ -547,7 +589,7 @@ detect_version <- function(cells) { file_desc <- sub("^EUCAST\\s+", "", version_text) file_desc <- trimws(sub(",?\\s*valid from.*$", "", file_desc, ignore.case = TRUE)) if (!nzchar(file_desc)) file_desc <- NA_character_ - + list(version = version, guideline = guideline, file_desc = file_desc) } @@ -556,15 +598,15 @@ detect_version <- function(cells) { # ============================================================================== parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { cells <- xlsx_cells(xlsx_path, sheets = sheet_name) - + # --- Detect version --- ver <- detect_version(cells) - + formats <- xlsx_formats(xlsx_path) - + # --- Parse rich text for cols A:I --- parsed <- parse_cells(cells |> filter(col >= 1, col <= 9), formats) - + # --- Statements extending MIC breakpoints to another species --- # The H. influenzae sheet states (v9.0-v16.1): "In the absence of specific # breakpoints, the H. influenzae MIC breakpoints can be applied to H. @@ -586,21 +628,25 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { } parse_log$mic_scope_extensions <- rbind( parse_log$mic_scope_extensions, - data.frame(guideline = ver$guideline, sheet = sheet_name, - from_text = paste(genus, m[3]), to_text = paste(genus, m[5]), statement = m[1]) + data.frame( + guideline = ver$guideline, sheet = sheet_name, + from_text = paste(genus, m[3]), to_text = paste(genus, m[5]), statement = m[1] + ) ) } else if (grepl("breakpoints can (also )?be applied to [A-Z]\\.", txt_squished)) { - stop("parse_sheet(): unrecognised statement about applying breakpoints to another species on sheet '", - sheet_name, "': ", substr(txt_squished, 1, 300)) + stop( + "parse_sheet(): unrecognised statement about applying breakpoints to another species on sheet '", + sheet_name, "': ", substr(txt_squished, 1, 300) + ) } } - + # --- Detect header rows (col 2 contains "MIC breakpoint") --- header_rows <- parsed |> filter(col == 2, grepl("MIC breakpoint", base_text, fixed = TRUE)) |> pull(row) |> sort() - + if (length(header_rows) == 0) { # Only acceptable when EUCAST states on the sheet itself that it has not # determined breakpoints, e.g. "EUCAST has not determined breakpoints for @@ -613,9 +659,11 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { if (length(no_bp_statement) > 0) { parse_log$no_breakpoints_determined <- rbind( parse_log$no_breakpoints_determined, - data.frame(guideline = ver$guideline, sheet = sheet_name, - organism = trimws(get_cell(parsed, 1, 1)$base), - statement = gsub("\\s+", " ", trimws(no_bp_statement[1]))) + data.frame( + guideline = ver$guideline, sheet = sheet_name, + organism = trimws(get_cell(parsed, 1, 1)$base), + statement = gsub("\\s+", " ", trimws(no_bp_statement[1])) + ) ) message(" EUCAST has not determined breakpoints in sheet '", sheet_name, "', skipping") return(NULL) @@ -631,10 +679,10 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { } stop("parse_sheet(): no 'MIC breakpoint' header row found in sheet '", sheet_name, "'") } - + # Sub-header rows: one row below each header (contains "S " pattern) sub_header_rows <- header_rows + 1 - + # --- Detect the notes column --- # On sheets that report both MIC and disk breakpoints, columns are laid out # as: A agent | B-D MIC (S/R/ATU) | E-H disk (dose/S/R/ATU) | I notes. @@ -648,20 +696,20 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { header_row_cols <- parsed |> filter(row %in% header_rows, col > 3, grepl("^Notes", base_text)) |> pull(col) - + notes_col <- if (length(header_row_cols) > 0) { # Normally consistent across header rows within a sheet; take the modal # (most frequent) value defensively in case of stray formatting cells. as.integer(names(sort(table(header_row_cols), decreasing = TRUE))[1]) } else { - 9L # fallback to the conventional position if no "Notes" label is found + 9L # fallback to the conventional position if no "Notes" label is found } - + # Data columns run up to (but excluding) the notes column; disk columns # are only present when notes_col > 5 (i.e. the sheet is wider than a # MIC-only layout of A:E). has_disk_columns <- notes_col > 5 - + # --- Detect sheet type: single-organism vs multi-organism --- # Multi-organism: col A at header rows contains "Antimicrobial agent". # This label alone is not sufficient, though: some genuinely single- @@ -675,33 +723,39 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { header_col_a <- parsed |> filter(row %in% header_rows, col == 1) |> pull(base_text) - + is_multi_organism <- length(header_rows) > 1 && all(grepl("Antimicrobial agent", header_col_a, fixed = TRUE)) - + # --- Detect organisms and their table ranges --- max_data_row <- max(parsed$row) - + if (is_multi_organism) { # Multi-organism: discover organism names above each header row. # Organism name is a standalone row (col A only, no data in B:H) # found between the previous table's end and this header row. - rows_with_b <- parsed |> filter(col == 2) |> pull(row) - + rows_with_b <- parsed |> + filter(col == 2) |> + pull(row) + tables <- list() for (i in seq_along(header_rows)) { hr <- header_rows[i] shr <- sub_header_rows[i] - + # Search backwards from the header row for the organism name search_from <- if (i == 1) 1 else sub_header_rows[i - 1] + 1 candidates <- parsed |> - filter(col == 1, row >= search_from, row < hr, - !row %in% rows_with_b) |> + filter( + col == 1, row >= search_from, row < hr, + !row %in% rows_with_b + ) |> # Exclude known non-organism text patterns - filter(!grepl("Antimicrobial|MIC determination|Disk diffusion|Breakpoints for|Expert Rules|For abbreviations|For species|Examples of|haze|Isolated|Ignore haemolysis|Numbered|Lettered|Medium:|Inoculum:|Incubation:|Reading:|Quality control:|See disk", - base_text)) - + filter(!grepl( + "Antimicrobial|MIC determination|Disk diffusion|Breakpoints for|Expert Rules|For abbreviations|For species|Examples of|haze|Isolated|Ignore haemolysis|Numbered|Lettered|Medium:|Inoculum:|Incubation:|Reading:|Quality control:|See disk", + base_text + )) + # Take the last candidate (closest to the header) if (nrow(candidates) > 0) { org_row <- candidates |> slice_max(row, n = 1) @@ -710,7 +764,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # Fallback: use sheet name organism <- sheet_name } - + # Data range: from sub_header + 1 to the row before the next organism # or end of data first_data <- shr + 1 @@ -731,15 +785,18 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { last_data <- last_data - 1 } } - + # Notes cell: notes_col, in the data range (usually at first_data, merged) - note_cell <- parsed |> filter(col == notes_col, row >= first_data, row <= last_data) |> + note_cell <- parsed |> + filter(col == notes_col, row >= first_data, row <= last_data) |> slice_min(row, n = 1) note_text <- if (nrow(note_cell) > 0) note_cell$full_text[1] else NA_character_ notes_parsed <- parse_notes_block(note_text) - - tables[[i]] <- list(organism = organism, first_data = first_data, - last_data = last_data, notes = notes_parsed) + + tables[[i]] <- list( + organism = organism, first_data = first_data, + last_data = last_data, notes = notes_parsed + ) } } else { # Single-organism: organism name from row 1, col A @@ -747,7 +804,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { organism <- if (nrow(org_cell) > 0) org_cell$base_text[1] else sheet_name # Clean: remove trailing * or whitespace organism <- gsub("[*]+$", "", trimws(organism)) - + # Notes: notes_col, within each class section # For single-organism sheets, there is one notes block per class. # Notes cell is at the header row or first data row of each class. @@ -756,7 +813,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { hr <- header_rows[i] shr <- sub_header_rows[i] first_data <- shr + 1 - + if (i < length(header_rows)) { last_data <- header_rows[i + 1] - 1 while (last_data >= first_data) { @@ -772,23 +829,27 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { last_data <- last_data - 1 } } - + # Notes: look in notes_col from header row through end of section note_cells <- parsed |> - filter(col == notes_col, row >= hr, row <= last_data, - !grepl("^Notes", base_text)) + filter( + col == notes_col, row >= hr, row <= last_data, + !grepl("^Notes", base_text) + ) note_text <- if (nrow(note_cells) > 0) { paste(note_cells$full_text, collapse = "\n") } else { NA_character_ } notes_parsed <- parse_notes_block(note_text) - - tables[[i]] <- list(organism = organism, first_data = first_data, - last_data = last_data, notes = notes_parsed) + + tables[[i]] <- list( + organism = organism, first_data = first_data, + last_data = last_data, notes = notes_parsed + ) } } - + # --- Extract breakpoint rows --- # Split agent name from any organism restriction in column A (see # split_agent_organism() for the formatting-based rationale). This is @@ -800,16 +861,16 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # "unexpected format" error at all. data_rows <- unique(unlist(map(tables, function(tbl) tbl$first_data:tbl$last_data))) agent_org_split <- split_agent_organism(cells, formats, sheet_name, valid_rows = data_rows) - + results <- vector("list", 500) idx <- 0L - + for (tbl in tables) { org <- tbl$organism - nl <- tbl$notes + nl <- tbl$notes prev_agent_base <- NA_character_ prev_disk_dose <- NA_character_ - + for (r in tbl$first_data:tbl$last_data) { agent <- get_cell(parsed, r, 1) # Skip category headers (rows that have col A text but no data). @@ -817,8 +878,10 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # columns; on MIC-only sheets col 5 is the notes column and must not # be treated as a data indicator. data_cols <- if (has_disk_columns) c(2, 3, 5, 6, 7) else c(2, 3) - has_any_data <- any(!is.na(vapply(data_cols, function(cc) get_cell(parsed, r, cc)$base, - character(1)))) + has_any_data <- any(!is.na(vapply( + data_cols, function(cc) get_cell(parsed, r, cc)$base, + character(1) + ))) if (is.na(agent$base)) { # A row with values but without an agent name would otherwise be # skipped silently, e.g. if a future table merged agent cells. Rows @@ -831,7 +894,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { next } if (!has_any_data) next - + # Agent name and any organism restriction encoded in its formatting # (see split_agent_organism()). Falls back to the raw cell text and # the table-level organism when this row has no split entry (e.g. @@ -843,25 +906,25 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { } else { org } - - mic_s <- get_cell(parsed, r, 2) - mic_r <- get_cell(parsed, r, 3) + + mic_s <- get_cell(parsed, r, 2) + mic_r <- get_cell(parsed, r, 3) mic_atu <- get_cell(parsed, r, 4) # Disk columns (E:H) only exist on sheets that report disk diffusion # breakpoints. On MIC-only sheets, col E is the notes column instead # of the disk dose, so disk cells must not be read from it there. if (has_disk_columns) { disk_dose_cell <- get_cell(parsed, r, 5) - disk_s <- get_cell(parsed, r, 6) - disk_r <- get_cell(parsed, r, 7) + disk_s <- get_cell(parsed, r, 6) + disk_r <- get_cell(parsed, r, 7) disk_atu <- get_cell(parsed, r, 8) } else { disk_dose_cell <- list(base = NA_character_, super = NA_character_) - disk_s <- list(base = NA_character_, super = NA_character_) - disk_r <- list(base = NA_character_, super = NA_character_) + disk_s <- list(base = NA_character_, super = NA_character_) + disk_r <- list(base = NA_character_, super = NA_character_) disk_atu <- list(base = NA_character_, super = NA_character_) } - + # EUCAST occasionally leaves the disk content empty in a row that # continues the previous row's agent, e.g. v16.1 "Amoxicillin- # clavulanic acid, F. nucleatum" directly below "Amoxicillin- @@ -883,7 +946,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { } prev_agent_base <- agent_base prev_disk_dose <- disk_dose - + # Notes must apply only to the method they were referenced from: MIC # superscripts resolve to the MIC row's note, disk superscripts to the # DISK row's note. Resolving both against the union of MIC and disk @@ -894,67 +957,71 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { combine_super <- function(...) { parts <- c(...) parts <- parts[!is.na(parts) & parts != ""] - if (length(parts) == 0) return(NA_character_) + if (length(parts) == 0) { + return(NA_character_) + } paste(parts, collapse = ",") } - mic_note <- resolve_notes(nl, combine_super(mic_s$super, agent$super), NA, abbreviations = abbreviations) + mic_note <- resolve_notes(nl, combine_super(mic_s$super, agent$super), NA, abbreviations = abbreviations) disk_note <- resolve_notes(nl, NA, combine_super(disk_s$super, agent$super), abbreviations = abbreviations) - + # MIC row s_val <- parse_bp(mic_s$base) r_val <- parse_bp(mic_r$base) if (!is.na(s_val) || !is.na(r_val)) { idx <- idx + 1L results[[idx]] <- tibble( - guideline = ver$guideline, - version = ver$version, - file_desc = ver$file_desc, - type = "human", - host = "human", - method = "MIC", - site = NA_character_, - mo = row_mo, + guideline = ver$guideline, + version = ver$version, + file_desc = ver$file_desc, + type = "human", + host = "human", + method = "MIC", + site = NA_character_, + mo = row_mo, table_organism = org, - rank_index = NA_integer_, - ab = agent_text, - disk_dose = NA_character_, + rank_index = NA_integer_, + ab = agent_text, + disk_dose = NA_character_, breakpoint_S = s_val, breakpoint_R = r_val, is_bracketed = is_bracketed_bp(mic_s$base) || is_bracketed_bp(mic_r$base), - note = mic_note + note = mic_note ) } - + # DISK row s_val <- parse_bp(disk_s$base) r_val <- parse_bp(disk_r$base) if (!is.na(s_val) || !is.na(r_val)) { idx <- idx + 1L results[[idx]] <- tibble( - guideline = ver$guideline, - version = ver$version, - file_desc = ver$file_desc, - type = "human", - host = "human", - method = "DISK", - site = NA_character_, - mo = row_mo, + guideline = ver$guideline, + version = ver$version, + file_desc = ver$file_desc, + type = "human", + host = "human", + method = "DISK", + site = NA_character_, + mo = row_mo, table_organism = org, - rank_index = NA_integer_, - ab = agent_text, - disk_dose = disk_dose, + rank_index = NA_integer_, + ab = agent_text, + disk_dose = disk_dose, breakpoint_S = s_val, breakpoint_R = r_val, is_bracketed = is_bracketed_bp(disk_s$base) || is_bracketed_bp(disk_r$base), - note = disk_note + note = disk_note ) } } } - - if (idx == 0) return(NULL) + + if (idx == 0) { + return(NULL) + } out <- bind_rows(results[seq_len(idx)]) - + # The non-species related PK-PD breakpoints (v9.0-v13.1) come with EUCAST's # condition for their use, e.g. "These breakpoints are used only when there # are no species-specific breakpoints or other recommendations (a dash or a @@ -962,7 +1029,7 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # the note of every PK-PD row. if (grepl("Non-species related", organism, ignore.case = TRUE)) { usage <- parsed$full_text[!is.na(parsed$full_text) & parsed$col == 1 & - grepl("^These breakpoints are used only when", trimws(parsed$full_text))] + grepl("^These breakpoints are used only when", trimws(parsed$full_text))] if (length(usage) != 1) { stop("parse_sheet(): the condition for using the PK-PD breakpoints was not found on sheet '", sheet_name, "'") } @@ -977,9 +1044,9 @@ parse_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # ============================================================================== parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { cells <- xlsx_cells(xlsx_path, sheets = "Topical agents") - + ver <- detect_version(cells) - + parsed <- parse_cells(cells, xlsx_formats(xlsx_path)) # v9.0 has a different, informational table here ("ECOFFs and systemic @@ -996,8 +1063,8 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { # on, disk contents in row 10, and organism rows (MIC row followed by a # zone diameter row) from row 11 on. Verify this rather than assume it. if (!identical(get_cell(parsed, 6, 1)$base, "Organisms") || - !identical(get_cell(parsed, 10, 2)$base, "Disk content") || - !grepl("^Screening cut-off values", coalesce(get_cell(parsed, 6, 2)$base, ""))) { + !identical(get_cell(parsed, 10, 2)$base, "Disk content") || + !grepl("^Screening cut-off values", coalesce(get_cell(parsed, 6, 2)$base, ""))) { stop("parse_topical_sheet(): unexpected layout of sheet 'Topical agents' in '", basename(xlsx_path), "'") } @@ -1011,15 +1078,18 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { ab_info[[as.character(col_idx)]] <- list( name = gsub("\r\n", " ", name_cell$base), name_super = name_cell$super, - dose = if (!is.na(dose_cell$base) && !dose_cell$base %in% c("-", "ND")) - paste0(dose_cell$base, " mcg") else NA_character_ + dose = if (!is.na(dose_cell$base) && !dose_cell$base %in% c("-", "ND")) { + paste0(dose_cell$base, " mcg") + } else { + NA_character_ + } ) } - + # Organism rows: col 1 from row 11 onwards, excluding "Notes" org_cells <- parsed |> filter(col == 1, row >= 11, !grepl("^Notes", base_text)) - + organism_pairs <- list() for (i in seq_len(nrow(org_cells))) { mic_row <- org_cells$row[i] @@ -1033,20 +1103,26 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { mic_row = mic_row, disk_row = disk_row ) } - + # Notes: find in col 1, rows after the last organism data row - last_data_row <- if (length(organism_pairs) > 0) - max(vapply(organism_pairs, function(x) x$disk_row, numeric(1))) else 0L + last_data_row <- if (length(organism_pairs) > 0) { + max(vapply(organism_pairs, function(x) x$disk_row, numeric(1))) + } else { + 0L + } notes_cells <- parsed |> filter(col == 1, row > last_data_row, !grepl("^Notes$", base_text)) - notes_text <- if (nrow(notes_cells) > 0) - paste(notes_cells$full_text, collapse = "\n") else "" + notes_text <- if (nrow(notes_cells) > 0) { + paste(notes_cells$full_text, collapse = "\n") + } else { + "" + } if (is.na(notes_text)) notes_text <- "" notes_lookup <- parse_notes_block(notes_text) - + results <- vector("list", 200) idx <- 0L - + for (od in organism_pairs) { for (col_idx in ab_cols) { cc <- as.character(col_idx) @@ -1055,7 +1131,7 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { ab_name <- ab_entry$name ab_super <- ab_entry$name_super ab_dose <- ab_entry$dose - + # MIC mic_cell <- get_cell(parsed, od$mic_row, col_idx) mic_val <- parse_bp(mic_cell$base) @@ -1063,7 +1139,7 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { all_supers <- c(mic_cell$super, ab_super) all_supers <- all_supers[!is.na(all_supers) & all_supers != ""] mic_note <- resolve_notes(notes_lookup, paste(all_supers, collapse = ","), NA, abbreviations = abbreviations) - + if (!is.na(mic_val)) { idx <- idx + 1L results[[idx]] <- tibble( @@ -1077,14 +1153,14 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { note = mic_note ) } - + # DISK disk_cell <- get_cell(parsed, od$disk_row, col_idx) disk_val <- parse_bp(disk_cell$base) all_supers_d <- c(disk_cell$super, ab_super) all_supers_d <- all_supers_d[!is.na(all_supers_d) & all_supers_d != ""] disk_note <- resolve_notes(notes_lookup, NA, paste(all_supers_d, collapse = ","), abbreviations = abbreviations) - + if (!is.na(disk_val)) { idx <- idx + 1L results[[idx]] <- tibble( @@ -1100,8 +1176,10 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { } } } - - if (idx == 0) return(NULL) + + if (idx == 0) { + return(NULL) + } bind_rows(results[seq_len(idx)]) } @@ -1122,11 +1200,11 @@ parse_topical_sheet <- function(xlsx_path, abbreviations = list()) { # diffusion method, so only method = "MIC" rows are produced. parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { cells <- xlsx_cells(xlsx_path, sheets = sheet_name) - + ver <- detect_version(cells) - + parsed <- parse_cells(cells, xlsx_formats(xlsx_path)) - + # --- Locate the "Antifungal agent" header row --- hdr_row <- parsed |> filter(col == 1, grepl("^Antifungal agent$", base_text)) |> @@ -1139,19 +1217,23 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { species_row <- hdr_row + 1L sub_row <- hdr_row + 2L first_data_row <- hdr_row + 3L - + # --- Detect species column-blocks from the sub-header row --- # Each block starts at a cell reading "S <=" (S \u2264) and runs through the # next "R >" and, where present, an immediately following "ATU" column, # stopping before the next "S <=" or the end of populated columns. - sub_cells <- parsed |> filter(row == sub_row) |> arrange(col) - s_cols <- sub_cells |> filter(grepl("^S\\b", base_text)) |> pull(col) - + sub_cells <- parsed |> + filter(row == sub_row) |> + arrange(col) + s_cols <- sub_cells |> + filter(grepl("^S\\b", base_text)) |> + pull(col) + if (length(s_cols) == 0) { message(" No species S/R columns found in sheet '", sheet_name, "', skipping") return(NULL) } - + species_blocks <- list() for (i in seq_along(s_cols)) { s_col <- s_cols[i] @@ -1160,7 +1242,7 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { # block is 3 columns wide rather than 2) atu_cell <- get_cell(parsed, sub_row, s_col + 2L) has_atu <- !is.na(atu_cell$base) && grepl("^ATU\\b", atu_cell$base) - + sp_cell <- get_cell(parsed, species_row, s_col) species <- if (!is.na(sp_cell$base)) gsub("\r?\n", " ", sp_cell$base) else NA_character_ if (is.na(species)) next @@ -1179,13 +1261,13 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { if (grepl("^[A-Z]\\.\\s", species) && substr(species, 1, 1) == substr(sheet_genus, 1, 1)) { species <- sub("^[A-Z]\\.", sheet_genus, species) } - + species_blocks[[length(species_blocks) + 1]] <- list( species = species, s_col = s_col, r_col = r_col, atu_col = if (has_atu) s_col + 2L else NA_integer_ ) } - + # --- Locate the end of the data block (row before the "Notes" label) --- notes_label_row <- parsed |> filter(col == 1, grepl("^Notes$", base_text), row > first_data_row) |> @@ -1195,7 +1277,7 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { } else { max(parsed$row) } - + # --- Notes block: one note per row below the "Notes" label, col A --- notes_text <- if (length(notes_label_row) > 0) { note_cells <- parsed |> @@ -1206,52 +1288,54 @@ parse_yeast_sheet <- function(xlsx_path, sheet_name, abbreviations = list()) { "" } notes_lookup <- parse_notes_block(notes_text) - + # --- Extract breakpoint rows: one per (agent, species) combination --- results <- vector("list", 1000) idx <- 0L - + for (r in first_data_row:last_data_row) { agent_cell <- get_cell(parsed, r, 1) if (is.na(agent_cell$base)) next agent_name <- gsub("\r?\n", " ", agent_cell$base) agent_super <- agent_cell$super - + for (sb in species_blocks) { s_cell <- get_cell(parsed, r, sb$s_col) r_cell <- get_cell(parsed, r, sb$r_col) - + s_val <- parse_bp(s_cell$base) r_val <- parse_bp(r_cell$base) if (is.na(s_val) && is.na(r_val)) next - + all_supers <- c(s_cell$super, r_cell$super, agent_super) all_supers <- all_supers[!is.na(all_supers) & all_supers != ""] note_text <- resolve_notes(notes_lookup, paste(all_supers, collapse = ","), NA, abbreviations = abbreviations) - + idx <- idx + 1L results[[idx]] <- tibble( - guideline = ver$guideline, - version = ver$version, - file_desc = ver$file_desc, - type = "human", - host = "human", - method = "MIC", - site = NA_character_, - mo = sb$species, + guideline = ver$guideline, + version = ver$version, + file_desc = ver$file_desc, + type = "human", + host = "human", + method = "MIC", + site = NA_character_, + mo = sb$species, table_organism = sb$species, - rank_index = NA_integer_, - ab = agent_name, - disk_dose = NA_character_, + rank_index = NA_integer_, + ab = agent_name, + disk_dose = NA_character_, breakpoint_S = s_val, breakpoint_R = r_val, is_bracketed = is_bracketed_bp(s_cell$base) || is_bracketed_bp(r_cell$base), - note = note_text + note = note_text ) } } - - if (idx == 0) return(NULL) + + if (idx == 0) { + return(NULL) + } bind_rows(results[seq_len(idx)]) } @@ -1289,19 +1373,21 @@ parse_workbook_notes <- function(xlsx_path) { all_sheets <- xlsx_sheet_names(xlsx_path) notes_sheet <- all_sheets[grepl("^\\s*(\\d+\\.\\s*)?Notes\\s*$", all_sheets, ignore.case = TRUE)] if (length(notes_sheet) != 1) { - stop("parse_workbook_notes(): expected exactly one 'Notes' sheet in '", basename(xlsx_path), - "', found: ", paste(notes_sheet, collapse = ", ")) + stop( + "parse_workbook_notes(): expected exactly one 'Notes' sheet in '", basename(xlsx_path), + "', found: ", paste(notes_sheet, collapse = ", ") + ) } txt <- parse_cells(xlsx_cells(xlsx_path, sheets = notes_sheet), xlsx_formats(xlsx_path)) |> filter(!is.na(full_text)) |> pull(full_text) |> gsub(pattern = "\u00A0", replacement = " ", fixed = TRUE) |> trimws() - + abbr_txt <- txt[grepl("^[A-Z]{2,5}\\s*=\\s*\\S", txt)] abbreviations <- as.list(trimws(sub("^[A-Z]{2,5}\\s*=\\s*", "", abbr_txt))) names(abbreviations) <- sub("\\s*=.*$", "", abbr_txt) - + bracket_definition <- NA_character_ def_brackets_sign <- txt[grepl("^\\(\\s*\\)\\s*=", txt)] def_own_note <- sub("^\\d+\\.\\s*", "", txt[grepl("^(\\d+\\.\\s*)?Breakpoints in brackets", txt)]) @@ -1314,7 +1400,7 @@ parse_workbook_notes <- function(xlsx_path) { bracket_definition <- sub("^.*?(Breakpoints in brackets)", "\\1", def_in_note[1]) } bracket_definition <- trimws(gsub("\\s+", " ", bracket_definition)) - + list(abbreviations = abbreviations, bracket_definition = bracket_definition) } @@ -1322,12 +1408,14 @@ parse_workbook_notes <- function(xlsx_path) { # 10. Parse all data sheets # ============================================================================== parse_workbook <- function(xlsx_path, - skip_sheets = c("Content", "Changes", "Notes", - "Guidance", "Dosages", - "Technical uncertainty", - "Non-species related breakpoints")) { + skip_sheets = c( + "Content", "Changes", "Notes", + "Guidance", "Dosages", + "Technical uncertainty", + "Non-species related breakpoints" + )) { all_sheets <- xlsx_sheet_names(xlsx_path) - + # Some workbooks (e.g. the antifungal AFST tables) number-prefix their # sheet names, e.g. "1. Notes", "7. Dosages", "6. Aspergillus ", and are # not always consistent in capitalisation or separators, e.g. @@ -1337,8 +1425,10 @@ parse_workbook <- function(xlsx_path, # (with underscores read as spaces) by whether a skip name occurs anywhere # in the (trimmed) sheet name instead. is_skipped <- vapply(gsub("_", " ", trimws(all_sheets)), function(s) { - any(vapply(skip_sheets, function(skip) grepl(tolower(skip), tolower(s), fixed = TRUE), - logical(1))) + any(vapply( + skip_sheets, function(skip) grepl(tolower(skip), tolower(s), fixed = TRUE), + logical(1) + )) }, logical(1)) # Hidden sheets are not part of the published tables: v14.0-v16.1 contain a # hidden, outdated copy of the v14.0 PK-PD table with editorial comments @@ -1348,9 +1438,9 @@ parse_workbook <- function(xlsx_path, message(" skipping hidden sheet(s): ", paste(hidden, collapse = ", ")) } data_sheets <- all_sheets[!is_skipped & !all_sheets %in% hidden] - + wb_notes <- parse_workbook_notes(xlsx_path) - + results <- list() for (s in data_sheets) { message("Parsing ", basename(xlsx_path), ": ", s) @@ -1378,21 +1468,24 @@ parse_workbook <- function(xlsx_path, message(" -> ", nrow(res), " rows") } } - + out <- bind_rows(results) - + # Add EUCAST's own definition of breakpoints in brackets to the note of # every bracketed row, so the caveat travels with the values in the raw # data (they are excluded from the clinical_breakpoints table). if (any(out$is_bracketed)) { if (is.na(wb_notes$bracket_definition)) { - stop("parse_workbook(): '", basename(xlsx_path), "' has bracketed breakpoints, ", - "but no definition of breakpoints in brackets was found on its Notes sheet") + stop( + "parse_workbook(): '", basename(xlsx_path), "' has bracketed breakpoints, ", + "but no definition of breakpoints in brackets was found on its Notes sheet" + ) } bracket_note <- paste0("[( )] ", wb_notes$bracket_definition) out$note[out$is_bracketed] <- ifelse(is.na(out$note[out$is_bracketed]), - bracket_note, - paste(bracket_note, out$note[out$is_bracketed], sep = " | ")) + bracket_note, + paste(bracket_note, out$note[out$is_bracketed], sep = " | ") + ) } out } @@ -1401,33 +1494,34 @@ parse_workbook <- function(xlsx_path, # 11. Run # ============================================================================== -breakpoint_files <- list.files(path = "data-raw", - pattern = "breakpoint.*table.*[.]xlsx$", - full.names = TRUE, - recursive = FALSE, - ignore.case = TRUE) +breakpoint_files <- list.files( + path = "data-raw", + pattern = "breakpoint.*table.*[.]xlsx$", + full.names = TRUE, + recursive = FALSE, + ignore.case = TRUE +) breakpoint_files <- breakpoint_files[breakpoint_files %unlike% "dosages"] breakpoints_eucast <- tibble() for (xlsx_path in breakpoint_files) { - message("Parsing ", basename(xlsx_path)) result <- suppressMessages(parse_workbook(xlsx_path)) - + # cat("\n============================\n") cat("Total rows:", nrow(result), "\n\n") # glimpse(result) - + # cat("\n--- Rows per sheet ---\n") # result |> count(sheet) |> print(n = 40) - + # cat("\n--- Rows per organism (top 20) ---\n") # result |> count(mo) |> arrange(desc(n)) |> head(20) |> print() - + breakpoints_eucast <- breakpoints_eucast |> bind_rows(result) - + # saveRDS(result, "/home/claude/eucast_all_breakpoints.rds") # write.csv(result, "/home/claude/eucast_all_breakpoints.csv", row.names = FALSE) # cat("\nSaved.\n") @@ -1436,8 +1530,10 @@ for (xlsx_path in breakpoint_files) { breakpoints_eucast_raw <- breakpoints_eucast if (!is.null(parse_log$disk_dose_inherited)) { - message("NOTE: ", nrow(parse_log$disk_dose_inherited), " disk breakpoint row(s) had no disk content and took it ", - "from the previous row of the same agent:") + message( + "NOTE: ", nrow(parse_log$disk_dose_inherited), " disk breakpoint row(s) had no disk content and took it ", + "from the previous row of the same agent:" + ) print(parse_log$disk_dose_inherited) } @@ -1505,10 +1601,10 @@ split_organism_from_ab <- function(ab, mo) { # Strip parenthetical content first to test for a bare comma. ab_no_parens <- gsub("\\([^)]*\\)", "", ab) has_bare_comma <- grepl(",", ab_no_parens) - + ab_clean <- ab mo_from_ab <- rep(NA_character_, length(ab)) - + for (i in which(has_bare_comma)) { # Split at the first bare comma (outside parentheses) # Simple approach: find the first comma in ab_no_parens, use its @@ -1516,23 +1612,24 @@ split_organism_from_ab <- function(ab, mo) { comma_pos <- regexpr(",", ab_no_parens[i]) if (comma_pos < 1) next before <- trimws(substr(ab[i], 1, comma_pos - 1)) - after <- trimws(substr(ab[i], comma_pos + 1, nchar(ab[i]))) + after <- trimws(substr(ab[i], comma_pos + 1, nchar(ab[i]))) if (!nzchar(after)) next - + # The part after the comma must look like an organism restriction, # not a drug-name continuation. Heuristic: starts with a genus-like # capitalised word (3+ letters), or starts with "other" / "except", # or contains known taxonomic suffixes (cocci, bacill, etc.) is_organism <- grepl( "^(other\\b|except\\b|[a-z]*cocc|[a-z]*bacill|[A-Z][a-z]{2,}\\b)", - after, ignore.case = FALSE, perl = TRUE + after, + ignore.case = FALSE, perl = TRUE ) if (!is_organism) next - + ab_clean[i] <- before mo_from_ab[i] <- after } - + tibble(ab_clean = ab_clean, mo_from_ab = mo_from_ab) } @@ -1637,8 +1734,10 @@ build_site <- function(sheet_site, route, qualifier) { q <- tolower(gsub("\\s+", " ", trimws(qualifier))) unknown <- unique(qualifier[!is.na(q) & !q %in% names(eucast_site_vocabulary)]) if (length(unknown) > 0) { - stop("build_site(): unknown indication qualifier(s), add to eucast_site_vocabulary: '", - paste(unknown, collapse = "', '"), "'") + stop( + "build_site(): unknown indication qualifier(s), add to eucast_site_vocabulary: '", + paste(unknown, collapse = "', '"), "'" + ) } indication <- unname(eucast_site_vocabulary[q]) vapply(seq_along(q), function(i) { @@ -1702,7 +1801,7 @@ expand_abbreviated_genus <- function(pieces, fallback_genus = NA_character_) { return(p) } if (!is.na(current_genus) && grepl("^[A-Z]\\.\\s", p) && - substr(p, 1, 1) == substr(current_genus, 1, 1)) { + substr(p, 1, 1) == substr(current_genus, 1, 1)) { return(sub("^[A-Z]\\.", current_genus, p)) } if (!is.na(current_genus) && grepl("^[a-z]+$", p)) { @@ -1728,7 +1827,7 @@ split_mo_list <- function(x, sheet_name = NA_character_) { # not matching the "Genus species" shape on at least one piece, and are # returned unchanged so as.mo() resolves them directly to their # species-group code. - + # Resolve each distinct sheet's own organism once (not per row -- x and # sheet_name are typically thousands of rows spanning ~40 distinct # sheets), to use as an abbreviation seed (see expand_abbreviated_genus()) @@ -1746,7 +1845,7 @@ split_mo_list <- function(x, sheet_name = NA_character_) { genus_lookup[sn] <- mo_genus(sheet_mo) } } - + map2(x, sheet_name, function(s, sn) { sheet_genus_fallback <- if (!is.na(sn) && sn %in% names(genus_lookup)) genus_lookup[[sn]] else NA_character_ split_mo_list_one(s, sheet_genus_fallback) @@ -1807,9 +1906,11 @@ split_mo_list_one <- function(s, sheet_genus_fallback) { pieces <- pieces[pieces != ""] pieces <- strip_mo_exceptions(pieces) pieces <- trimws(gsub("\\s*\\(\\s*\\)\\s*$", "", pieces)) # drop now-empty "()" leftovers - - if (length(pieces) <= 1) return(s_clean) - + + if (length(pieces) <= 1) { + return(s_clean) + } + # Only treat as a genuine multi-taxon list if every piece looks like a # taxon on its own: either an abbreviated genus initial ("E.", "P." # followed by a species epithet), a capitalised word of 2+ letters @@ -1844,8 +1945,10 @@ split_mo_list_one <- function(s, sheet_genus_fallback) { # whole cell falls back to the original, unsplit string. looks_like_taxon <- grepl("^[A-Z]\\.\\s", pieces) | grepl("^[A-Z][a-z]+(\\.|\\s|$)", pieces) | grepl("^[a-z]+$", pieces) | grepl("^coagulase-negative staphylococci$", pieces, ignore.case = TRUE) - if (!all(looks_like_taxon)) return(s_clean) - + if (!all(looks_like_taxon)) { + return(s_clean) + } + # Reject a piece containing more than one capitalised, 3+ letter word # (a rough genus/proper-name detector): this signals two taxa were # fused into one piece because a separating comma was missing in the @@ -1856,10 +1959,14 @@ split_mo_list_one <- function(s, sheet_genus_fallback) { # as.mo("Citrobacter spp. Klebsiella spp.") resolves to Citrobacter # enshiensis) -- so this must be caught before as.mo() ever sees it, # by falling back to the unsplit original string. - genus_like_count <- vapply(pieces, function(p) length(gregexpr("\\b[A-Z][a-z]{2,}\\b", p, perl = TRUE)[[1]]), - integer(1)) - if (any(genus_like_count > 1)) return(s_clean) - + genus_like_count <- vapply( + pieces, function(p) length(gregexpr("\\b[A-Z][a-z]{2,}\\b", p, perl = TRUE)[[1]]), + integer(1) + ) + if (any(genus_like_count > 1)) { + return(s_clean) + } + expand_abbreviated_genus(pieces, fallback_genus = sheet_genus_fallback) } @@ -1908,10 +2015,14 @@ breakpoints_eucast <- breakpoints_eucast_raw |> is_blocker = is_bracketed | (!is_bp_value(breakpoint_S) & !is_bp_value(breakpoint_R)), note = case_when( is_bracketed ~ blocking_note(paste0("EUCAST gives breakpoints in brackets only (S ", breakpoint_S, ", R ", breakpoint_R, ")"), note), - is_blocker ~ blocking_note(paste0("EUCAST gives ", - coalesce(bp_placeholder_meaning[toupper(trimws(breakpoint_S))], - bp_placeholder_meaning[toupper(trimws(breakpoint_R))], - "no value")), note), + is_blocker ~ blocking_note(paste0( + "EUCAST gives ", + coalesce( + bp_placeholder_meaning[toupper(trimws(breakpoint_S))], + bp_placeholder_meaning[toupper(trimws(breakpoint_R))], + "no value" + ) + ), note), TRUE ~ note ), # placeholders become NA; a row keeps a breakpoint as long as one of S and @@ -1929,8 +2040,10 @@ breakpoints_eucast <- breakpoints_eucast_raw |> ab_mo_fix <- split_organism_from_ab(breakpoints_eucast$ab, breakpoints_eucast$mo) fixed_rows <- !is.na(ab_mo_fix$mo_from_ab) if (any(fixed_rows)) { - message("NOTE: ", sum(fixed_rows), " row(s) had an organism restriction embedded in the ab column; ", - "splitting into ab + mo (e.g. '", breakpoints_eucast$ab[which(fixed_rows)[1]], "').") + message( + "NOTE: ", sum(fixed_rows), " row(s) had an organism restriction embedded in the ab column; ", + "splitting into ab + mo (e.g. '", breakpoints_eucast$ab[which(fixed_rows)[1]], "')." + ) breakpoints_eucast$ab[fixed_rows] <- ab_mo_fix$ab_clean[fixed_rows] breakpoints_eucast$mo[fixed_rows] <- ab_mo_fix$mo_from_ab[fixed_rows] } @@ -1939,8 +2052,10 @@ if (any(fixed_rows)) { breakpoints_eucast <- breakpoints_eucast |> mutate(note = case_when( trimws(gsub("\\s+", " ", mo)) != "MRSA" ~ note, - is_blocker ~ trimws(paste0(trimws(sub("^(\\[No breakpoint\\][^|]*).*$", "\\1", note)), " | [Scope] EUCAST gives this for MRSA ", - sub("^\\[No breakpoint\\][^|]*", "", note))), + is_blocker ~ trimws(paste0( + trimws(sub("^(\\[No breakpoint\\][^|]*).*$", "\\1", note)), " | [Scope] EUCAST gives this for MRSA ", + sub("^\\[No breakpoint\\][^|]*", "", note) + )), is.na(note) ~ "[Scope] EUCAST gives this breakpoint for MRSA", TRUE ~ paste("[Scope] EUCAST gives this breakpoint for MRSA", note, sep = " | ") )) @@ -1981,7 +2096,8 @@ extract_mo_exceptions <- function(x) { pieces <- pieces[nzchar(pieces)] if (nzchar(genus)) { pieces <- ifelse(grepl("^[A-Z]\\.\\s*[a-z]", pieces) & substr(pieces, 1, 1) == substr(genus, 1, 1), - sub("^[A-Z]\\.\\s*", paste0(genus, " "), pieces), pieces) + sub("^[A-Z]\\.\\s*", paste0(genus, " "), pieces), pieces + ) } pieces })) @@ -1994,9 +2110,11 @@ exceptions <- breakpoints_eucast |> tidyr::unnest(exception) if (nrow(exceptions) > 0) { exception_blockers <- exceptions |> - mutate(note = blocking_note(paste0("EUCAST excludes ", exception, " ('", gsub("\\s+", " ", mo), "')"), NA_character_), - mo = exception, site = NA_character_, is_blocker = TRUE, - breakpoint_S = NA_character_, breakpoint_R = NA_character_) |> + mutate( + note = blocking_note(paste0("EUCAST excludes ", exception, " ('", gsub("\\s+", " ", mo), "')"), NA_character_), + mo = exception, site = NA_character_, is_blocker = TRUE, + breakpoint_S = NA_character_, breakpoint_R = NA_character_ + ) |> select(-exception) breakpoints_eucast <- bind_rows(breakpoints_eucast, exception_blockers) } @@ -2042,21 +2160,25 @@ mo_uncertainties_eucast <- AMR:::AMR_env$mo_uncertainties # does not cover the group, as happened with "S. anginosus group". unmapped_groups <- breakpoints_eucast |> distinct(mo_text, mo) |> - filter(grepl("group|anaerob|MRSA|MSSA|-negative|-positive|cocci\\b|bacteria\\b", mo_text, ignore.case = TRUE), - !grepl("spp\\.$", mo_text), - # a text that is exactly the name of the taxon it resolved to, such - # as "Peptostreptococcus anaerobius", is not a descriptive group - tolower(mo_text) != tolower(mo_name(mo, keep_synonyms = TRUE)), - !mo_rank(mo, keep_synonyms = TRUE) %in% "species group", - !mo %in% unlist(eucast_organism_groups)) + filter( + grepl("group|anaerob|MRSA|MSSA|-negative|-positive|cocci\\b|bacteria\\b", mo_text, ignore.case = TRUE), + !grepl("spp\\.$", mo_text), + # a text that is exactly the name of the taxon it resolved to, such + # as "Peptostreptococcus anaerobius", is not a descriptive group + tolower(mo_text) != tolower(mo_name(mo, keep_synonyms = TRUE)), + !mo_rank(mo, keep_synonyms = TRUE) %in% "species group", + !mo %in% unlist(eucast_organism_groups) + ) if (nrow(unmapped_groups) > 0) { print(unmapped_groups |> mutate(name = mo_name(mo, keep_synonyms = TRUE))) stop(nrow(unmapped_groups), " descriptive organism group(s) resolved to a single taxon, add them to eucast_organism_groups") } if (nrow(mo_uncertainties_eucast) > 0) { - message("NOTE: ", nrow(mo_uncertainties_eucast), " organism string(s) were resolved with uncertainty ", - "(matched below the default confidence, or against multiple candidates). ", - "Inspect `mo_uncertainties_eucast` and cross-check against the source sheet before trusting these rows.") + message( + "NOTE: ", nrow(mo_uncertainties_eucast), " organism string(s) were resolved with uncertainty ", + "(matched below the default confidence, or against multiple candidates). ", + "Inspect `mo_uncertainties_eucast` and cross-check against the source sheet before trusting these rows." + ) } breakpoints_eucast <- breakpoints_eucast |> @@ -2100,9 +2222,11 @@ ab_previously_coerced_eucast <- AMR:::AMR_env$ab_previously_coerced |> filter(tolower(trimws(x_bak)) %in% tolower(trimws(unique(breakpoints_eucast$ab_input)))) breakpoints_eucast <- breakpoints_eucast |> select(-ab_input) if (nrow(ab_previously_coerced_eucast) > 0) { - message("NOTE: ", nrow(ab_previously_coerced_eucast), " antibiotic string(s) were resolved via approximate/previous-coercion ", - "matching rather than an exact code or name match. ", - "Inspect `ab_previously_coerced_eucast` and cross-check against the source sheet before trusting these rows.") + message( + "NOTE: ", nrow(ab_previously_coerced_eucast), " antibiotic string(s) were resolved via approximate/previous-coercion ", + "matching rather than an exact code or name match. ", + "Inspect `ab_previously_coerced_eucast` and cross-check against the source sheet before trusting these rows." + ) } # Restrictions in the notes ---- @@ -2185,8 +2309,10 @@ taxon_contains <- function(broad, narrow) { broad %in% AMR::microorganisms.groups$mo_group[AMR::microorganisms.groups$mo == as.character(narrow)] } to_blocker <- function(df, reason) { - df |> mutate(is_blocker = TRUE, breakpoint_S = NA_character_, breakpoint_R = NA_character_, - site = NA_character_, note = blocking_note(reason, note)) + df |> mutate( + is_blocker = TRUE, breakpoint_S = NA_character_, breakpoint_R = NA_character_, + site = NA_character_, note = blocking_note(reason, note) + ) } for (i in seq_len(nrow(eucast_note_rules))) { rule <- eucast_note_rules[i, ] @@ -2211,9 +2337,11 @@ for (i in seq_len(nrow(eucast_note_rules))) { breakpoints_eucast[rows, ] <- to_blocker(breakpoints_eucast[rows, ], reason_of(rows)) } else if (rule$action == "site") { breakpoints_eucast$site[governed] <- if_else(is.na(breakpoints_eucast$site[governed]), rule$target, - if_else(grepl(rule$target, breakpoints_eucast$site[governed], fixed = TRUE), - breakpoints_eucast$site[governed], - paste(breakpoints_eucast$site[governed], rule$target, sep = ", "))) + if_else(grepl(rule$target, breakpoints_eucast$site[governed], fixed = TRUE), + breakpoints_eucast$site[governed], + paste(breakpoints_eucast$site[governed], rule$target, sep = ", ") + ) + ) } else if (rule$action == "restrict") { rows <- which(governed & !breakpoints_eucast$is_blocker) for (r in rows) { @@ -2269,8 +2397,10 @@ derived_rows <- list() # (only the rows of the test itself, e.g. "Gentamicin (test for high-level # aminoglycoside resistance)"; other aminoglycosides refer to the same note) hlar_rows <- breakpoints_eucast |> - filter(grepl("Negative test: Isolates with (gentamicin|streptomycin) MIC \u2264", note), - ab %in% c("GEN", "STR1"), site %in% "Screen") + filter( + grepl("Negative test: Isolates with (gentamicin|streptomycin) MIC \u2264", note), + ab %in% c("GEN", "STR1"), site %in% "Screen" + ) for (r in seq_len(nrow(hlar_rows))) { row <- hlar_rows[r, ] # the note of a row can also hold the test for the other agent, so the @@ -2285,38 +2415,46 @@ for (r in seq_len(nrow(hlar_rows))) { value <- if (row$method == "MIC") neg[2] else neg[3] if (!nzchar(value)) next # e.g. no zone diameter given in this version derived_rows[[length(derived_rows) + 1]] <- row |> - mutate(ab = as.ab(if (agent == "gentamicin") "GEH" else "STH"), site = NA_character_, is_blocker = FALSE, - breakpoint_S = value, breakpoint_R = value, - note = sub("^\\[No breakpoint\\][^|]*(\\| )?", "", note)) + mutate( + ab = as.ab(if (agent == "gentamicin") "GEH" else "STH"), site = NA_character_, is_blocker = FALSE, + breakpoint_S = value, breakpoint_R = value, + note = sub("^\\[No breakpoint\\][^|]*(\\| )?", "", note) + ) } fox_mic <- breakpoints_eucast |> - filter(ab == "FOX", method == "MIC", - grepl("S\\. aureus and S\\. lugdunensis with cefoxitin MIC values >[0-9]+ mg/L and S\\. saprophyticus with cefoxitin MIC values >[0-9]+ mg/L are methicillin resistant", note)) |> + filter( + ab == "FOX", method == "MIC", + grepl("S\\. aureus and S\\. lugdunensis with cefoxitin MIC values >[0-9]+ mg/L and S\\. saprophyticus with cefoxitin MIC values >[0-9]+ mg/L are methicillin resistant", note) + ) |> distinct(guideline, .keep_all = TRUE) for (r in seq_len(nrow(fox_mic))) { row <- fox_mic[r, ] m <- regmatches(row$note, regexec("S\\. aureus and S\\. lugdunensis with cefoxitin MIC values >([0-9]+) mg/L and S\\. saprophyticus with cefoxitin MIC values >([0-9]+) mg/L", row$note))[[1]] for (sp in list(c("B_STPHY_AURS", m[2]), c("B_STPHY_LGDN", m[2]), c("B_STPHY_SPRP", m[3]))) { derived_rows[[length(derived_rows) + 1]] <- row |> - mutate(mo = as.mo(sp[1]), site = "Screen", is_blocker = FALSE, breakpoint_S = sp[2], breakpoint_R = sp[2], - note = sub("^\\[No breakpoint\\][^|]*(\\| )?", "", note)) + mutate( + mo = as.mo(sp[1]), site = "Screen", is_blocker = FALSE, breakpoint_S = sp[2], breakpoint_R = sp[2], + note = sub("^\\[No breakpoint\\][^|]*(\\| )?", "", note) + ) } } if (length(derived_rows) > 0) { derived <- bind_rows(derived_rows) - message("NOTE: ", nrow(derived), " screening breakpoint row(s) taken from the notes (", - paste(unique(paste(derived$guideline, derived$ab, derived$method)), collapse = ", "), ")") + message( + "NOTE: ", nrow(derived), " screening breakpoint row(s) taken from the notes (", + paste(unique(paste(derived$guideline, derived$ab, derived$method)), collapse = ", "), ")" + ) breakpoints_eucast <- bind_rows(breakpoints_eucast, derived) } cons_rule <- "If coagulase-negative staphylococci are not identified to species level, use zone diameter breakpoints S\u2265([0-9]+), R<([0-9]+) mm" for (g in unique(breakpoints_eucast$guideline[grepl(cons_rule, breakpoints_eucast$note) & breakpoints_eucast$ab == "FOX" & breakpoints_eucast$method == "DISK"])) { rows <- which(breakpoints_eucast$guideline == g & breakpoints_eucast$ab == "FOX" & breakpoints_eucast$method == "DISK" & - breakpoints_eucast$mo == "B_STPHY_CONS" & !breakpoints_eucast$is_blocker) + breakpoints_eucast$mo == "B_STPHY_CONS" & !breakpoints_eucast$is_blocker) if (length(rows) != 1) stop("expected one cefoxitin zone diameter row for coagulase-negative staphylococci in ", g) m <- regmatches(breakpoints_eucast$note[rows], regexec(cons_rule, breakpoints_eucast$note[rows]))[[1]] if (length(m) != 3 || m[2] != m[3]) stop("unexpected cefoxitin rule for unidentified CoNS in ", g) own_rows <- unique(as.character(breakpoints_eucast$mo[breakpoints_eucast$guideline == g & breakpoints_eucast$ab == "FOX" & - breakpoints_eucast$method == "DISK" & breakpoints_eucast$mo != "B_STPHY_CONS"])) + breakpoints_eucast$method == "DISK" & breakpoints_eucast$mo != "B_STPHY_CONS"])) members <- setdiff(as.character(AMR::microorganisms.groups$mo[AMR::microorganisms.groups$mo_group == "B_STPHY_CONS"]), own_rows) species_rows <- breakpoints_eucast[rep(rows, length(members)), ] species_rows$mo <- as.mo(members) @@ -2363,8 +2501,10 @@ breakpoints_eucast <- bind_rows(breakpoints_eucast, scope_blockers) if (!is.null(parse_log$mic_scope_extensions)) { scope_ext <- parse_log$mic_scope_extensions |> distinct() |> - mutate(from_mo = as.character(as.mo(from_text, info = FALSE)), - to_mo = as.character(as.mo(to_text, info = FALSE))) + mutate( + from_mo = as.character(as.mo(from_text, info = FALSE)), + to_mo = as.character(as.mo(to_text, info = FALSE)) + ) if (anyNA(scope_ext$from_mo) || anyNA(scope_ext$to_mo) || any(c(scope_ext$from_mo, scope_ext$to_mo) == "UNKNOWN")) { print(scope_ext) stop("species in a MIC scope statement could not be resolved, see above") @@ -2373,18 +2513,25 @@ if (!is.null(parse_log$mic_scope_extensions)) { filter(method == "MIC") |> mutate(mo_chr = as.character(mo)) |> inner_join(scope_ext |> select(guideline, sheet, mo_chr = from_mo, to_mo, statement), - by = c("guideline", "sheet", "mo_chr")) |> - mutate(mo = as.mo(to_mo), - note = case_when( - is_blocker ~ trimws(paste0(trimws(sub("^(\\[No breakpoint\\][^|]*).*$", "\\1", note)), " | [Scope] ", statement, " ", - sub("^\\[No breakpoint\\][^|]*", "", note))), - is.na(note) ~ paste0("[Scope] ", statement), - TRUE ~ paste0("[Scope] ", statement, " | ", note) - )) |> + by = c("guideline", "sheet", "mo_chr") + ) |> + mutate( + mo = as.mo(to_mo), + note = case_when( + is_blocker ~ trimws(paste0( + trimws(sub("^(\\[No breakpoint\\][^|]*).*$", "\\1", note)), " | [Scope] ", statement, " ", + sub("^\\[No breakpoint\\][^|]*", "", note) + )), + is.na(note) ~ paste0("[Scope] ", statement), + TRUE ~ paste0("[Scope] ", statement, " | ", note) + ) + ) |> select(-mo_chr, -to_mo, -statement) |> anti_join(breakpoints_eucast |> filter(method == "MIC"), by = c("guideline", "mo", "ab", "method", "site")) - message("NOTE: ", nrow(scope_rows), " MIC breakpoint row(s) added following EUCAST's statement(s): ", - paste(unique(paste0(scope_ext$guideline, " ", scope_ext$from_text, " -> ", scope_ext$to_text)), collapse = ", ")) + message( + "NOTE: ", nrow(scope_rows), " MIC breakpoint row(s) added following EUCAST's statement(s): ", + paste(unique(paste0(scope_ext$guideline, " ", scope_ext$from_text, " -> ", scope_ext$to_text)), collapse = ", ") + ) breakpoints_eucast <- bind_rows(breakpoints_eucast, scope_rows) } @@ -2458,10 +2605,10 @@ breakpoints_eucast <- breakpoints_eucast |> })) |> # Greek symbols and EM dash symbols are not allowed by CRAN, so replace them with ASCII: mutate(disk_dose = disk_dose %>% - gsub("\u03bc", "mc", ., fixed = TRUE) %>% # this is 'mu', \u03bc - gsub("\u00b5", "mc", ., fixed = TRUE) %>% # this is 'micro', \u00b5 (yes, they look the same) - gsub("\u2013", "-", ., fixed = TRUE) %>% - gsub("(?<=\\d)(?=[a-zA-Z])", " ", ., perl = TRUE)) # keep a space after a number, e.g. "1mcg" to "1 mcg" + gsub("\u03bc", "mc", ., fixed = TRUE) %>% # this is 'mu', \u03bc + gsub("\u00b5", "mc", ., fixed = TRUE) %>% # this is 'micro', \u00b5 (yes, they look the same) + gsub("\u2013", "-", ., fixed = TRUE) %>% + gsub("(?<=\\d)(?=[a-zA-Z])", " ", ., perl = TRUE)) # keep a space after a number, e.g. "1mcg" to "1 mcg" # Every row must have resolved to an organism and an antimicrobial. as.mo() # returns the code "UNKNOWN" rather than NA when it finds no match at all @@ -2470,8 +2617,12 @@ breakpoints_eucast <- breakpoints_eucast |> # the parsing above (e.g. a new free-text pattern), never a genuine ambiguity. # The only exception is the non-species related PK-PD table, which is coded # as "UNKNOWN" on purpose (as.sir() uses it as the last resort). -unresolved_mo <- breakpoints_eucast |> filter(is.na(mo) | (mo == "UNKNOWN" & !grepl("^PK", sheet))) |> distinct(sheet, ref_tbl, ab, mo_text) -unresolved_ab <- breakpoints_eucast |> filter(is.na(ab)) |> distinct(sheet, ref_tbl) +unresolved_mo <- breakpoints_eucast |> + filter(is.na(mo) | (mo == "UNKNOWN" & !grepl("^PK", sheet))) |> + distinct(sheet, ref_tbl, ab, mo_text) +unresolved_ab <- breakpoints_eucast |> + filter(is.na(ab)) |> + distinct(sheet, ref_tbl) if (nrow(unresolved_mo) > 0 || nrow(unresolved_ab) > 0) { print(unresolved_mo, n = Inf) print(unresolved_ab, n = Inf) @@ -2544,17 +2695,21 @@ if (nrow(missing_dose) > 0) { # S breakpoint, and only screening tests may lack an R breakpoint. bad_rows <- breakpoints_eucast |> filter((is_blocker & (!is.na(breakpoint_S) | !is.na(breakpoint_R) | !grepl("^\\[No breakpoint\\] ", note))) | - (!is_blocker & (is.na(breakpoint_S) | (is.na(breakpoint_R) & !grepl("screen", site, ignore.case = TRUE))))) + (!is_blocker & (is.na(breakpoint_S) | (is.na(breakpoint_R) & !grepl("screen", site, ignore.case = TRUE))))) if (nrow(bad_rows) > 0) { print(bad_rows |> select(guideline, mo, ab, method, site, breakpoint_S, breakpoint_R, ref_tbl, note), n = Inf) stop(nrow(bad_rows), " rows are neither a complete breakpoint nor a proper blocking row, see above") } breakpoints_eucast <- breakpoints_eucast |> - select(guideline, type, host, method, site, mo, rank_index, ab, ref_tbl, - disk_dose, breakpoint_S, breakpoint_R, uti, is_SDD, note) + select( + guideline, type, host, method, site, mo, rank_index, ab, ref_tbl, + disk_dose, breakpoint_S, breakpoint_R, uti, is_SDD, note + ) -breakpoints_eucast |> count(guideline) |> print() +breakpoints_eucast |> + count(guideline) |> + print() # Which guideline years these workbooks cover, separately for bacteria and # fungi: the bacterial tables and the antifungal (AFST) tables are separate @@ -2577,8 +2732,9 @@ merge_eucast_breakpoints <- function(clinical_breakpoints, breakpoints_eucast, e # the antifungal tables scope_of <- function(mo, ab) { if_else(mo_kingdom(mo, keep_synonyms = TRUE) %in% "Fungi" | - (as.character(mo) == "UNKNOWN" & ab_group(ab) %in% "Antifungals"), - "Fungi", "Bacteria") + (as.character(mo) == "UNKNOWN" & ab_group(ab) %in% "Antifungals"), + "Fungi", "Bacteria" + ) } # the organisms in the EUCAST rows must agree with the type of workbook # they come from @@ -2595,9 +2751,11 @@ merge_eucast_breakpoints <- function(clinical_breakpoints, breakpoints_eucast, e filter(type == "human") kept <- clinical_breakpoints |> anti_join(replaced, by = names(clinical_breakpoints)) - message("Replacing ", nrow(replaced), " WHONET-based human EUCAST rows (", - paste(sort(unique(paste(replaced$guideline, replaced$scope))), collapse = ", "), - ") with ", nrow(breakpoints_eucast), " rows from the EUCAST workbooks") + message( + "Replacing ", nrow(replaced), " WHONET-based human EUCAST rows (", + paste(sort(unique(paste(replaced$guideline, replaced$scope))), collapse = ", "), + ") with ", nrow(breakpoints_eucast), " rows from the EUCAST workbooks" + ) # rows without any breakpoint are only kept as the blocking rows of the # EUCAST workbooks (see "Blocking rows") kept |> diff --git a/data-raw/_reproduction_scripts/reproduction_of_intrinsic_resistant.R b/data-raw/_reproduction_scripts/reproduction_of_intrinsic_resistant.R index e84f5f0e3..9681cbb40 100644 --- a/data-raw/_reproduction_scripts/reproduction_of_intrinsic_resistant.R +++ b/data-raw/_reproduction_scripts/reproduction_of_intrinsic_resistant.R @@ -38,12 +38,12 @@ for (i in seq_len(nrow(antimicrobials))) { } int_resis <- interpretive_rules(int_resis, - interpretive_rules_df = subset( - AMR:::INTERPRETIVE_RULES_DF, - reference.rule_group %like% "Expected.* Phenotypes" - ), - overwrite = TRUE, - info = FALSE + interpretive_rules_df = subset( + AMR:::INTERPRETIVE_RULES_DF, + reference.rule_group %like% "Expected.* Phenotypes" + ), + overwrite = TRUE, + info = FALSE ) int_resis2 <- int_resis[, sapply(int_resis, function(x) any(!is.sir(x) | x == "R")), drop = FALSE] %>% @@ -63,8 +63,10 @@ all(int_resis3$ab %in% antimicrobials$ab) intrinsic_resistant <- int_resis3 -intrinsic_resistant$version <- paste(EUCAST_VERSION_EXPECTED_PHENOTYPES$`1.2`$title, - EUCAST_VERSION_EXPECTED_PHENOTYPES$`1.2`$version_txt) +intrinsic_resistant$version <- paste( + EUCAST_VERSION_EXPECTED_PHENOTYPES$`1.2`$title, + EUCAST_VERSION_EXPECTED_PHENOTYPES$`1.2`$version_txt +) usethis::use_data(intrinsic_resistant, internal = FALSE, overwrite = TRUE, version = 2, compress = "xz") rm(intrinsic_resistant) diff --git a/data-raw/_reproduction_scripts/reproduction_of_microorganisms.R b/data-raw/_reproduction_scripts/reproduction_of_microorganisms.R index 3f1f4ed3f..083fe4dca 100644 --- a/data-raw/_reproduction_scripts/reproduction_of_microorganisms.R +++ b/data-raw/_reproduction_scripts/reproduction_of_microorganisms.R @@ -53,7 +53,7 @@ # 5. Go to BacDive data base for the oxygen tolerance and cell shape. # Go to https://bacdive.dsmz.de/advsearch, filter 'Oxygen tolerance' or # 'Cell shape' on "*" and click Submit and click on the 'Download tabel as CSV' button. -# Last time these worked: +# Last time these worked: # - wget -O bacdive_oxygen_tolerance.csv 'https://bacdive.dsmz.de/advsearch/csv?fg%5B1%5D%5Bgc%5D=OR&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfd%5D=Oxygen+tolerance&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfv%5D=%2A&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfvd%5D=oxygen_tolerance-oxygen_tol-4' # - wget -O bacdive_cell_shape.csv 'https://bacdive.dsmz.de/advsearch/csv?fg%5B1%5D%5Bgc%5D=OR&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfd%5D=Cell+shape&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfv%5D=%2A&fg%5B1%5D%5Bfl%5D%5B1%5D%5Bfvd%5D=cell_morphology-cell_shape-2' # 6. Set these locations to the paths where the files are (or set the environment variable @@ -231,12 +231,12 @@ get_author_year <- function(ref) { # normalise Unicode to NFC (precomposed) so iconv can transliterate authors2 <- stringi::stri_trans_nfc(ref) # fix known encoding errors in source data - authors2 <- gsub("\u0092", "'", authors2) # Windows-1252 right quote - authors2 <- gsub("\u0091", "'", authors2) # Windows-1252 left quote - authors2 <- gsub("\u04AB", "\u00E7", authors2) # Cyrillic cedilla to Latin cedilla - authors2 <- gsub("[\u0400-\u04FF]", "", authors2) # remove remaining Cyrillic characters + authors2 <- gsub("\u0092", "'", authors2) # Windows-1252 right quote + authors2 <- gsub("\u0091", "'", authors2) # Windows-1252 left quote + authors2 <- gsub("\u04AB", "\u00E7", authors2) # Cyrillic cedilla to Latin cedilla + authors2 <- gsub("[\u0400-\u04FF]", "", authors2) # remove remaining Cyrillic characters authors2 <- iconv(authors2, from = "UTF-8", to = "ASCII//TRANSLIT") - + authors2 <- gsub(" ?\\(Approved Lists [0-9]+\\) ?", " ", authors2) authors2 <- gsub(" +", " ", authors2) authors2 <- trimws(authors2) @@ -251,15 +251,15 @@ get_author_year <- function(ref) { # replace parentheses with emend. to get the latest authors # authors2 <- gsub("(", " emend. ", authors2, fixed = TRUE) # authors2 <- gsub(")", "", authors2, fixed = TRUE) - + # remove any remaining parentheses authors2 <- gsub("[()]", "", authors2) authors2 <- gsub(" +", " ", authors2) authors2 <- trimws(authors2) - + # strip emend. and everything after it to retain the combination authority authors2 <- gsub(" ?emend[.]?.*", "", authors2) - + # get year from last 4 digits lastyear <- as.integer(gsub(".*([0-9]{4})$", "\\1", authors2)) # can never be later than now @@ -359,16 +359,19 @@ get_lpsn_and_author <- function(rank, name, tries = 3) { for (try in seq_len(tries)) { # be polite to the LPSN server Sys.sleep(0.1) - page_txt <- tryCatch({ - resp <- curl::curl_fetch_memory(url) - if (resp$status_code == 200) { - read_html(rawToChar(resp$content)) - } else { - # 404 is a valid answer, other codes (e.g. 429, 503) are worth a retry - not_found <- resp$status_code == 404 - NULL - } - }, error = function(e) NULL) + page_txt <- tryCatch( + { + resp <- curl::curl_fetch_memory(url) + if (resp$status_code == 200) { + read_html(rawToChar(resp$content)) + } else { + # 404 is a valid answer, other codes (e.g. 429, 503) are worth a retry + not_found <- resp$status_code == 404 + NULL + } + }, + error = function(e) NULL + ) if (!is.null(page_txt) || not_found) { break } @@ -1431,35 +1434,43 @@ taxonomy_gbif <- taxonomy_gbif %>% mutate( name_parts = strsplit(fullname, " "), rank_num = match(rank, rank_levels, nomatch = NA_integer_), - + # Step 1: parse the synonym's own name components from fullname kingdom = if_else(is.na(kingdom) & rank == "kingdom", - sapply(name_parts, `[`, 1L), kingdom), + sapply(name_parts, `[`, 1L), kingdom + ), phylum = if_else(is.na(phylum) & rank == "phylum", - sapply(name_parts, `[`, 1L), phylum), + sapply(name_parts, `[`, 1L), phylum + ), class = if_else(is.na(class) & rank == "class", - sapply(name_parts, `[`, 1L), class), + sapply(name_parts, `[`, 1L), class + ), order = if_else(is.na(order) & rank == "order", - sapply(name_parts, `[`, 1L), order), + sapply(name_parts, `[`, 1L), order + ), family = if_else(is.na(family) & rank == "family", - sapply(name_parts, `[`, 1L), family), + sapply(name_parts, `[`, 1L), family + ), genus = if_else(is.na(genus) & rank %in% c("genus", "species", "subspecies"), - sapply(name_parts, `[`, 1L), genus), + sapply(name_parts, `[`, 1L), genus + ), species = if_else(is.na(species) & rank %in% c("species", "subspecies"), - sapply(name_parts, `[`, 2L), species), + sapply(name_parts, `[`, 2L), species + ), subspecies = if_else(is.na(subspecies) & rank == "subspecies", - sapply(name_parts, `[`, 3L), subspecies), - + sapply(name_parts, `[`, 3L), subspecies + ), + # Step 2: fill domain and kingdom from the accepted name, only for ranks ABOVE the synonym's rank - domain = if_else(is.na(domain) & !is.na(acc_domain), acc_domain, domain), + domain = if_else(is.na(domain) & !is.na(acc_domain), acc_domain, domain), kingdom = if_else(is.na(kingdom) & !is.na(acc_kingdom) & rank_num > 1L, acc_kingdom, kingdom), # Step 3: phylum to family may only come from the accepted name if that is in the same genus # (or if the synonym is a genus or higher), otherwise from the synonym's own genus (see Step 4) same_lineage = rank_num <= 6L | (!is.na(acc_genus) & genus == acc_genus), - phylum = if_else(is.na(phylum) & !is.na(acc_phylum) & rank_num > 2L & same_lineage, acc_phylum, phylum), - class = if_else(is.na(class) & !is.na(acc_class) & rank_num > 3L & same_lineage, acc_class, class), - order = if_else(is.na(order) & !is.na(acc_order) & rank_num > 4L & same_lineage, acc_order, order), - family = if_else(is.na(family) & !is.na(acc_family) & rank_num > 5L & same_lineage, acc_family, family) + phylum = if_else(is.na(phylum) & !is.na(acc_phylum) & rank_num > 2L & same_lineage, acc_phylum, phylum), + class = if_else(is.na(class) & !is.na(acc_class) & rank_num > 3L & same_lineage, acc_class, class), + order = if_else(is.na(order) & !is.na(acc_order) & rank_num > 4L & same_lineage, acc_order, order), + family = if_else(is.na(family) & !is.na(acc_family) & rank_num > 5L & same_lineage, acc_family, family) ) %>% # Step 4: species in another genus than their accepted name get the taxonomy of their own genus left_join(accepted_genera, by = c("domain", "genus")) %>% @@ -1603,7 +1614,7 @@ taxonomy <- taxonomy_lpsn %>% ungroup() %>% # only keep unique species per domain distinct(domain, fullname, .keep_all = TRUE) %>% - arrange(fullname) %>% + arrange(fullname) %>% select(fullname, everything()) # get missing entries from existing microorganisms data set @@ -1845,12 +1856,12 @@ taxonomy <- taxonomy %>% filter(.resolve != "remove") %>% mutate( is_genus_rank = genus != "" & !is.na(best_domain) & rank %in% c("genus", "species", "subspecies"), - domain = if_else(genus != "" & !is.na(best_domain), best_domain, domain), + domain = if_else(genus != "" & !is.na(best_domain), best_domain, domain), kingdom = if_else(genus != "" & !is.na(best_kingdom) & best_kingdom != "", best_kingdom, kingdom), - phylum = if_else(is_genus_rank & best_phylum != "", best_phylum, phylum), - class = if_else(is_genus_rank & best_class != "", best_class, class), - order = if_else(is_genus_rank & best_order != "", best_order, order), - family = if_else(is_genus_rank & best_family != "", best_family, family) + phylum = if_else(is_genus_rank & best_phylum != "", best_phylum, phylum), + class = if_else(is_genus_rank & best_class != "", best_class, class), + order = if_else(is_genus_rank & best_order != "", best_order, order), + family = if_else(is_genus_rank & best_family != "", best_family, family) ) %>% select(-starts_with("best_"), -is_genus_rank, -.resolve) rm(genus_candidates, best_domain) @@ -1932,7 +1943,8 @@ taxonomy <- taxonomy %>% rank == "kingdom" ~ kingdom, rank == "domain" ~ domain, TRUE ~ paste(genus, species, subspecies) # already trimmed 7 lines up - ))) %>% + )) + ) %>% # keep the record of the most authoritative source: LPSN > MycoBank > GBIF > other mutate(sprio = source_prio(source)) %>% arrange(fullname, sprio) %>% @@ -2049,7 +2061,7 @@ fill_from_parent <- function(data, taxonomy, child_rank, parent_rank) { arrange(.data[[child_rank]], sprio) %>% distinct(.data[[child_rank]], .keep_all = TRUE) %>% select(all_of(c(child_rank, parent_rank))) - + data %>% rows_update(lookup, by = child_rank, unmatched = "ignore") } @@ -2091,7 +2103,7 @@ for (rank_name in c("family", "order", "class", "phylum", "kingdom")) { filter(rank != rank_name, .data[[rank_name]] != "") %>% distinct(.data[[rank_name]]) %>% pull() - + n_before <- nrow(taxonomy) taxonomy <- taxonomy %>% filter(!(rank == rank_name & !fullname %in% has_children)) @@ -2112,7 +2124,7 @@ gbif_bacteria <- which( taxonomy$domain == "Bacteria" & taxonomy$source %in% c("GBIF", "manually added", "inferred") & (taxonomy$rank %in% c("phylum", "class", "order", "family", "genus") | - (taxonomy$rank == "species" & taxonomy$genus %in% lpsn_genera)) + (taxonomy$rank == "species" & taxonomy$genus %in% lpsn_genera)) ) message(length(gbif_bacteria)) added <- 0 @@ -2175,18 +2187,19 @@ taxonomy$domain[taxonomy$domain == "(unknown kingdom)"] <- "(unknown domain)" # e.g., phylum above class "Bacilli" was still "Firmicutes" in 2023, should be "Bacillota" per LPSN rank_pairs <- tibble( child = c("genus", "family", "order", "class", "phylum"), - parent = c("family", "order", "class", "phylum", "kingdom") + parent = c("family", "order", "class", "phylum", "kingdom") ) for (d in unique(taxonomy$domain[taxonomy$domain != "(unknown domain)"])) { src <- if (d == "Fungi") "MycoBank" else "LPSN" message("Fixing GBIF taxonomy for domain ", d, " based on ", src, "...", - appendLF = FALSE) - + appendLF = FALSE + ) + for (i in seq_len(nrow(rank_pairs))) { - child_col <- rank_pairs$child[i] + child_col <- rank_pairs$child[i] parent_col <- rank_pairs$parent[i] - + # build lookup: for each child value in this domain, get the parent from the authoritative source lookup <- taxonomy %>% filter( @@ -2194,13 +2207,13 @@ for (d in unique(taxonomy$domain[taxonomy$domain != "(unknown domain)"])) { source == src, .data[[child_col]] != "" ) - + # special case: kingdom lookup excludes "Bacteria"/"Archaea" as kingdom values if (parent_col == "kingdom") { lookup <- lookup %>% filter(!kingdom %in% c("Bacteria", "Archaea")) } - + lookup <- lookup %>% # an empty parent in the authoritative source must never erase a known parent from GBIF # (e.g. MycoBank has no family for Mucor, GBIF has Mucoraceae) @@ -2209,7 +2222,7 @@ for (d in unique(taxonomy$domain[taxonomy$domain != "(unknown domain)"])) { arrange(rank != child_col) %>% distinct(.data[[child_col]], .keep_all = TRUE) %>% select(all_of(c(child_col, parent_col))) - + # join and overwrite taxonomy <- taxonomy %>% rows_update( @@ -2236,13 +2249,13 @@ taxonomy <- taxonomy %>% rank = case_when( rank %in% c("(unknown rank)", "species group") ~ rank, subspecies != "" ~ "subspecies", - species != "" ~ "species", - genus != "" ~ "genus", - family != "" ~ "family", - order != "" ~ "order", - class != "" ~ "class", - phylum != "" ~ "phylum", - kingdom != "" ~ "kingdom", + species != "" ~ "species", + genus != "" ~ "genus", + family != "" ~ "family", + order != "" ~ "order", + class != "" ~ "class", + phylum != "" ~ "phylum", + kingdom != "" ~ "kingdom", .default = NA_character_ ) ) @@ -2292,10 +2305,10 @@ taxonomy <- taxonomy %>% ) ) taxonomy$rank[taxonomy$kingdom == taxonomy$domain & - taxonomy$rank == "kingdom" & - taxonomy$domain %in% c("Bacteria", "Archaea")] <- "domain" + taxonomy$rank == "kingdom" & + taxonomy$domain %in% c("Bacteria", "Archaea")] <- "domain" taxonomy$kingdom[taxonomy$kingdom == taxonomy$domain & - taxonomy$rank == "domain"] <- "" + taxonomy$rank == "domain"] <- "" taxonomy <- taxonomy %>% arrange(fullname, source_prio(source), lpsn, mycobank, desc(gbif)) %>% @@ -2306,33 +2319,33 @@ taxonomy <- taxonomy %>% taxonomy <- taxonomy %>% mutate( .parent_fullname = case_when( - rank == "subspecies" ~ paste(genus, species), - rank == "species" ~ genus, - rank == "genus" & family != "" ~ family, - rank == "genus" & order != "" ~ order, - rank == "genus" & class != "" ~ class, - rank == "genus" & phylum != "" ~ phylum, - rank == "genus" ~ domain, - rank == "family" & order != "" ~ order, - rank == "family" & class != "" ~ class, - rank == "family" & phylum != "" ~ phylum, - rank == "family" ~ domain, - rank == "order" & class != "" ~ class, - rank == "order" & phylum != "" ~ phylum, - rank == "order" ~ domain, - rank == "class" & phylum != "" ~ phylum, - rank == "class" ~ domain, - rank == "phylum" & kingdom != "" ~ kingdom, - rank == "phylum" ~ domain, - rank == "kingdom" ~ domain, + rank == "subspecies" ~ paste(genus, species), + rank == "species" ~ genus, + rank == "genus" & family != "" ~ family, + rank == "genus" & order != "" ~ order, + rank == "genus" & class != "" ~ class, + rank == "genus" & phylum != "" ~ phylum, + rank == "genus" ~ domain, + rank == "family" & order != "" ~ order, + rank == "family" & class != "" ~ class, + rank == "family" & phylum != "" ~ phylum, + rank == "family" ~ domain, + rank == "order" & class != "" ~ class, + rank == "order" & phylum != "" ~ phylum, + rank == "order" ~ domain, + rank == "class" & phylum != "" ~ phylum, + rank == "class" ~ domain, + rank == "phylum" & kingdom != "" ~ kingdom, + rank == "phylum" ~ domain, + rank == "kingdom" ~ domain, .default = NA_character_ ) ) # Build lookup vectors: fullname -> id -lpsn_lookup <- setNames(taxonomy$lpsn, taxonomy$fullname) +lpsn_lookup <- setNames(taxonomy$lpsn, taxonomy$fullname) mycobank_lookup <- setNames(taxonomy$mycobank, taxonomy$fullname) -gbif_lookup <- setNames(taxonomy$gbif, taxonomy$fullname) +gbif_lookup <- setNames(taxonomy$gbif, taxonomy$fullname) taxonomy <- taxonomy %>% mutate( @@ -2676,12 +2689,12 @@ genus_ref <- taxonomy %>% filter(genus != "") %>% group_by(genus) %>% summarise( - .ref_domain = coalesce(first(domain[rank == "genus"]), first(domain)), - .ref_kingdom = coalesce(first(kingdom[rank == "genus"]), first(kingdom)), - .ref_phylum = consensus_value(phylum, rank, source, "genus"), - .ref_class = consensus_value(class, rank, source, "genus"), - .ref_order = consensus_value(order, rank, source, "genus"), - .ref_family = consensus_value(family, rank, source, "genus"), + .ref_domain = coalesce(first(domain[rank == "genus"]), first(domain)), + .ref_kingdom = coalesce(first(kingdom[rank == "genus"]), first(kingdom)), + .ref_phylum = consensus_value(phylum, rank, source, "genus"), + .ref_class = consensus_value(class, rank, source, "genus"), + .ref_order = consensus_value(order, rank, source, "genus"), + .ref_family = consensus_value(family, rank, source, "genus"), .groups = "drop" ) @@ -2689,13 +2702,13 @@ taxonomy <- taxonomy %>% left_join(genus_ref, by = "genus") %>% mutate( across_flag = rank != "genus" & genus != "" & fullname %unlike% "unknown", - domain = if_else(across_flag, .ref_domain, domain), - kingdom = if_else(across_flag, .ref_kingdom, kingdom), + domain = if_else(across_flag, .ref_domain, domain), + kingdom = if_else(across_flag, .ref_kingdom, kingdom), # also fill the genus record itself where it has no value - phylum = if_else((across_flag | (rank == "genus" & phylum == "")) & .ref_phylum != "", .ref_phylum, phylum), - class = if_else((across_flag | (rank == "genus" & class == "")) & .ref_class != "", .ref_class, class), - order = if_else((across_flag | (rank == "genus" & order == "")) & .ref_order != "", .ref_order, order), - family = if_else((across_flag | (rank == "genus" & family == "")) & .ref_family != "", .ref_family, family) + phylum = if_else((across_flag | (rank == "genus" & phylum == "")) & .ref_phylum != "", .ref_phylum, phylum), + class = if_else((across_flag | (rank == "genus" & class == "")) & .ref_class != "", .ref_class, class), + order = if_else((across_flag | (rank == "genus" & order == "")) & .ref_order != "", .ref_order, order), + family = if_else((across_flag | (rank == "genus" & family == "")) & .ref_family != "", .ref_family, family) ) %>% select(-starts_with(".ref_"), -across_flag) %>% arrange(fullname) @@ -2809,10 +2822,14 @@ plasmodia <- data.frame( status = "accepted", prevalence = 1.25 ) %>% - mutate(rank = ifelse(species == "", "genus", - ifelse(subspecies != "", "subspecies", - "species")), - fullname = trimws(paste(genus, species, subspecies))) + mutate( + rank = ifelse(species == "", "genus", + ifelse(subspecies != "", "subspecies", + "species" + ) + ), + fullname = trimws(paste(genus, species, subspecies)) + ) taxonomy <- taxonomy %>% filter(genus != "Plasmodium") %>% @@ -2898,8 +2915,10 @@ if (nrow(isospora_belli) == 1 && !"Cystoisospora belli" %in% taxonomy$fullname) status = "synonym", source = "manually added", ref = "", - across(any_of(c("lpsn", "lpsn_parent", "lpsn_renamed_to", "mycobank", "mycobank_parent", "mycobank_renamed_to", - "gbif", "gbif_parent")), ~NA_character_), + across(any_of(c( + "lpsn", "lpsn_parent", "lpsn_renamed_to", "mycobank", "mycobank_parent", "mycobank_renamed_to", + "gbif", "gbif_parent" + )), ~NA_character_), gbif_renamed_to = isospora_belli$gbif ) ) %>% @@ -3186,9 +3205,11 @@ taxonomy <- taxonomy %>% distinct(fullname, .keep_all = TRUE) %>% compute_prevalence() message("Restored ", nrow(restored), " released taxa (", nrow(taxonomy) - n_before, " records including their parents)") -rm(restore_registry, restore_renames, renamed, missing_released, homonyms_removed, retirements, new_retirements, source_names, older_homonyms, older_homonym_retirements, - released_records, genus_domains, conflicting, conflict_retirements, kingdom_of_phylum, lineage_of_genus, restored, - rank_order) +rm( + restore_registry, restore_renames, renamed, missing_released, homonyms_removed, retirements, new_retirements, source_names, older_homonyms, older_homonym_retirements, + released_records, genus_domains, conflicting, conflict_retirements, kingdom_of_phylum, lineage_of_genus, restored, + rank_order +) # Link orthographic variants without a current name ---------------------------------------------- @@ -3495,8 +3516,10 @@ if (nrow(new_retirements) > 0) { na = "" ) } -rm(no_current_name, lpsn_outcome, prokaryotes, to_remove, has_kept_children, outcome, removed, registry_now, - retirements_now, new_retirements) +rm( + no_current_name, lpsn_outcome, prokaryotes, to_remove, has_kept_children, outcome, removed, registry_now, + retirements_now, new_retirements +) # Fix genera that are synonyms while they contain accepted species -------------------------------- @@ -3617,8 +3640,10 @@ for (r in c("family", "order", "class", "phylum")) { new_record$status <- "accepted" new_record$source <- "manually added" new_record$ref <- NA_character_ - for (col in intersect(c("lpsn", "lpsn_parent", "lpsn_renamed_to", "mycobank", "mycobank_parent", - "mycobank_renamed_to", "gbif", "gbif_parent", "gbif_renamed_to"), colnames(new_record))) { + for (col in intersect(c( + "lpsn", "lpsn_parent", "lpsn_renamed_to", "mycobank", "mycobank_parent", + "mycobank_renamed_to", "gbif", "gbif_parent", "gbif_renamed_to" + ), colnames(new_record))) { new_record[[col]] <- NA_character_ } for (lf in c(r, lower)) new_record[[lf]] <- "" @@ -4022,7 +4047,7 @@ taxonomy <- taxonomy %>% filter(mo %unlike% "__") # keep unique taxonomy taxonomy <- bind_rows( - taxonomy %>% + taxonomy %>% filter(fullname %like% "unknown"), taxonomy %>% distinct(across(domain:subspecies), .keep_all = TRUE) %>% @@ -4489,7 +4514,8 @@ bacdive_shape_other <- taxonomy %>% morphology = case_when( fullname %like% "coccus" ~ "cocci", fullname %in% taxonomy$fullname[taxonomy$order %in% c("Enterobacterales", "Caryophanales", "Lactobacillales")] ~ morphology, - TRUE ~ paste("likely", morphology)) + TRUE ~ paste("likely", morphology) + ) ) %>% select(fullname, morphology) %>% distinct(fullname, .keep_all = TRUE) @@ -4883,7 +4909,7 @@ if (!is.null(mo_integrity_report(integrity_issues))) { # sets refer only to current names, except intrinsic_resistant, which lists synonyms on purpose) fix_old_mos <- function(dataset, new_ref, drop = FALSE, col = "mo", current = TRUE) { before <- dataset - + # the data sets contain codes of the development version mo_names <- microorganisms_dev$fullname[match(before[[col]], microorganisms_dev$mo)] # names of the development version that do not exist anymore, with their name in the new data set @@ -4922,20 +4948,24 @@ save_df <- function(...) { ) } microorganisms.codes <- fix_old_mos(AMR::microorganisms.codes, taxonomy, drop = TRUE) -if (!identical(microorganisms.codes, AMR::microorganisms.codes)) save_df(microorganisms.codes); rm(microorganisms.codes) +if (!identical(microorganisms.codes, AMR::microorganisms.codes)) save_df(microorganisms.codes) +rm(microorganisms.codes) clinical_breakpoints <- fix_old_mos(AMR::clinical_breakpoints, taxonomy) -if (!identical(clinical_breakpoints, AMR::clinical_breakpoints)) save_df(clinical_breakpoints); rm(clinical_breakpoints) +if (!identical(clinical_breakpoints, AMR::clinical_breakpoints)) save_df(clinical_breakpoints) +rm(clinical_breakpoints) example_isolates <- fix_old_mos(AMR::example_isolates, taxonomy) -if (!identical(example_isolates, AMR::example_isolates)) save_df(example_isolates); rm(example_isolates) +if (!identical(example_isolates, AMR::example_isolates)) save_df(example_isolates) +rm(example_isolates) # (intrinsic_resistant is not updated here: it is computed from the new data set afterwards, see # data-raw/_reproduction_scripts/reproduction_of_intrinsic_resistant.R) microorganisms.groups <- fix_old_mos(AMR::microorganisms.groups, taxonomy, col = "mo") microorganisms.groups <- fix_old_mos(microorganisms.groups, taxonomy, col = "mo_group") -if (!identical(microorganisms.groups, AMR::microorganisms.groups)) save_df(microorganisms.groups); rm(microorganisms.groups) +if (!identical(microorganisms.groups, AMR::microorganisms.groups)) save_df(microorganisms.groups) +rm(microorganisms.groups) # *** Save to package *** ------------------------------------------------------------------------- diff --git a/data-raw/_reproduction_scripts/run_after_microorganisms_build.R b/data-raw/_reproduction_scripts/run_after_microorganisms_build.R index 813e2b456..1f67765ac 100644 --- a/data-raw/_reproduction_scripts/run_after_microorganisms_build.R +++ b/data-raw/_reproduction_scripts/run_after_microorganisms_build.R @@ -56,12 +56,12 @@ step("microorganisms.groups", "data-raw/_reproduction_scripts/reproduction_of_mi # 2. EUCAST breakpoints, merged into clinical_breakpoints (as at the end of reproduction_of_clinical_breakpoints.R) merge_script <- tempfile(fileext = ".R") writeLines(c( - 'devtools::load_all()', - 'library(dplyr, warn.conflicts = FALSE)', + "devtools::load_all()", + "library(dplyr, warn.conflicts = FALSE)", 'source("data-raw/_reproduction_scripts/reproduction_of_clinical_breakpoints_eucast.R")', - '# the rank of an organism can have changed with the new data set', - 'breakpoints_now <- clinical_breakpoints %>%', - ' mutate(rank_index = case_when(', + "# the rank of an organism can have changed with the new data set", + "breakpoints_now <- clinical_breakpoints %>%", + " mutate(rank_index = case_when(", ' mo_rank(mo, keep_synonyms = TRUE) %like% "(infra|sub)" ~ 1,', ' mo_rank(mo, keep_synonyms = TRUE) == "species" ~ 2,', ' mo_rank(mo, keep_synonyms = TRUE) == "species group" ~ 2.5,', @@ -69,9 +69,9 @@ writeLines(c( ' mo_rank(mo, keep_synonyms = TRUE) == "family" ~ 4,', ' mo_rank(mo, keep_synonyms = TRUE) == "order" ~ 5,', ' mo != "UNKNOWN" ~ 6,', - ' TRUE ~ 7', + " TRUE ~ 7", ' ) + if_else(is.na(breakpoint_S) & is.na(breakpoint_R) & guideline %like% "EUCAST", 0.1, 0))', - 'clinical_breakpoints <- merge_eucast_breakpoints(breakpoints_now, breakpoints_eucast, eucast_coverage)', + "clinical_breakpoints <- merge_eucast_breakpoints(breakpoints_now, breakpoints_eucast, eucast_coverage)", 'usethis::use_data(clinical_breakpoints, overwrite = TRUE, compress = "xz", version = 2)' ), merge_script) step("clinical_breakpoints (EUCAST)", merge_script) diff --git a/data-raw/wisca_test/parameters_amr.R b/data-raw/wisca_test/parameters_amr.R index 770e6cf80..aab80b06f 100644 --- a/data-raw/wisca_test/parameters_amr.R +++ b/data-raw/wisca_test/parameters_amr.R @@ -1,622 +1,622 @@ # Copyright (c) [2022] [Larisse Bolton] -wisca_params <- function(x, antibiotic_in, pathogen_in, analysis, exclude, isolate_first, susceptible_I, infection_in,infection_full){ +wisca_params <- function(x, antibiotic_in, pathogen_in, analysis, exclude, isolate_first, susceptible_I, infection_in, infection_full) { # Isolate the type of syndrome you would like to investigate - if (infection_in != ""){ - x <- x %>% - filter(infection_type == infection_in) + if (infection_in != "") { + x <- x %>% + filter(infection_type == infection_in) } else { x <- x } - + # Define antibiotic regimens under investigation from "antimicrobials" dataset - antibiotic_df <- function(abo_in){ - mono_ab <- str_to_title(abo_in[!str_detect(abo_in,fixed("+"))]) #identify single antimicrobials - + antibiotic_df <- function(abo_in) { + mono_ab <- str_to_title(abo_in[!str_detect(abo_in, fixed("+"))]) # identify single antimicrobials + mono_antibiotic_x_set <- data.frame() - for (j in 1:length(mono_ab)){ #generate dataframe for single antimicrobials with abbreviations and fullnames + for (j in 1:length(mono_ab)) { # generate dataframe for single antimicrobials with abbreviations and fullnames mono_x_name <- ab_name(mono_ab[j], only_first = TRUE) mono_x_ab <- as.ab(mono_x_name) - mono_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% mono_x_ab)[,c("ab","name")] - mono_antibiotic_x_set <- rbind.data.frame(mono_antibiotic_x_set,mono_antibiotic, row.names = NULL) + mono_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% mono_x_ab)[, c("ab", "name")] + mono_antibiotic_x_set <- rbind.data.frame(mono_antibiotic_x_set, mono_antibiotic, row.names = NULL) } - names(mono_antibiotic_x_set) <- c("ab1","namerx") + names(mono_antibiotic_x_set) <- c("ab1", "namerx") mono_antibiotic_x_set$ab1 <- as.character(mono_antibiotic_x_set$ab1) - - - comb_ab <- str_to_title(abo_in[str_detect(abo_in,fixed("+"))]) #identify combination regimens - if(length(comb_ab) > 0){ - - comb_antibiotic_x_total <- data.frame() - for (jj in 1:length(comb_ab)){ #generate dataframe for combination antimicrobials with abbreviations and fullnames + + + comb_ab <- str_to_title(abo_in[str_detect(abo_in, fixed("+"))]) # identify combination regimens + if (length(comb_ab) > 0) { + comb_antibiotic_x_total <- data.frame() + for (jj in 1:length(comb_ab)) { # generate dataframe for combination antimicrobials with abbreviations and fullnames comb_ab_sep <- unlist(str_split(comb_ab[jj], fixed("+"))) comb_antibiotic_x_pre <- data.frame(x = rep(0, times = 2), y = rep(0, times = 2)) - for (ii in 1:length(comb_ab_sep)){ + for (ii in 1:length(comb_ab_sep)) { comb_x_name <- ab_name(comb_ab_sep[ii], only_first = TRUE) comb_x_ab <- as.ab(comb_x_name) - comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% comb_x_ab)[,c("ab","name")] - comb_antibiotic_x_pre[ii,] <- comb_antibiotic + comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% comb_x_ab)[, c("ab", "name")] + comb_antibiotic_x_pre[ii, ] <- comb_antibiotic names(comb_antibiotic_x_pre) <- names(comb_antibiotic) } comb_antibiotic_x <- pivot_wider(comb_antibiotic_x_pre, names_from = name, values_from = ab) comb_antibiotic_x_set <- comb_antibiotic_x %>% mutate(namerx = str_flatten(str_c(names(comb_antibiotic_x), collapse = "+"))) - new_ab <- sapply(X = as.character(1:(ncol(comb_antibiotic_x_set)-1)), FUN = function(x){str_flatten(str_c("ab",x))}) - names(comb_antibiotic_x_set)[1:(ncol(comb_antibiotic_x_set)-1)] <- new_ab - names(comb_antibiotic_x_set)[ncol(comb_antibiotic_x_set)]<- "namerx" - comb_antibiotic_x_total <-rbind.data.frame(comb_antibiotic_x_total, comb_antibiotic_x_set, row.names = NULL) + new_ab <- sapply(X = as.character(1:(ncol(comb_antibiotic_x_set) - 1)), FUN = function(x) { + str_flatten(str_c("ab", x)) + }) + names(comb_antibiotic_x_set)[1:(ncol(comb_antibiotic_x_set) - 1)] <- new_ab + names(comb_antibiotic_x_set)[ncol(comb_antibiotic_x_set)] <- "namerx" + comb_antibiotic_x_total <- rbind.data.frame(comb_antibiotic_x_total, comb_antibiotic_x_set, row.names = NULL) } - - - - antibiotic_x_set <- full_join(x= mono_antibiotic_x_set,y = comb_antibiotic_x_total) #generate full antibiotic dataset for analysis - antibiotic_x_set <- select(antibiotic_x_set, starts_with("ab"),starts_with("namerx")) - + + + antibiotic_x_set <- full_join(x = mono_antibiotic_x_set, y = comb_antibiotic_x_total) # generate full antibiotic dataset for analysis + antibiotic_x_set <- select(antibiotic_x_set, starts_with("ab"), starts_with("namerx")) } else { antibiotic_x_set <- mono_antibiotic_x_set } return(antibiotic_x_set) } - - antibiotic_rx <- antibiotic_df(abo_in = antibiotic_in) #full dataframe with antibiotic regimens - - - #if the dataset names are in different format to built-in set - if ((length(names(x)[sapply(names(x), function(x1){any(str_detect(x1,"[A-Z]$") & !(x1 %in% antimicrobials$ab) & !(x1 %in% antimicrobials$name))})])>0)){ - names(x)[which(str_detect(names(x),"[a-z]+") & !str_detect(names(x),"fullname") & !str_detect(names(x),"mo"))] <- str_to_upper(names(x)[which(str_detect(names(x),"[a-z]+") & !str_detect(names(x),"fullname") & !str_detect(names(x),"mo"))]) - names(x)[which(str_detect(names(x),"^ORG") & str_detect(names(x),"NAME$"))] <- "fullname" - names(x)[which(str_detect(names(x),"^ORG") & str_detect(names(x),"CODE$"))] <- "organism" - names(x)[which(str_detect(names(x),"^DATE") & str_detect(names(x),"CULTURE$") | str_detect(names(x),"^DATE") & str_detect(names(x),"SPECIMEN"))] <- "date" - names(x)[which(str_detect(names(x),"^EPISODE") | str_detect(names(x),"^PAT"))] <- "patient" - names(x)[which(str_detect(names(x),"^WARD"))] <- "ward" - names(x)[which(str_detect(names(x),"^HOSP"))] <- "hospital" - names(x)[which(str_detect(names(x),"DATE") & (str_detect(names(x),"REGISTRATION")))] <- "date" + antibiotic_rx <- antibiotic_df(abo_in = antibiotic_in) # full dataframe with antibiotic regimens + + + # if the dataset names are in different format to built-in set + if ((length(names(x)[sapply(names(x), function(x1) { + any(str_detect(x1, "[A-Z]$") & !(x1 %in% antimicrobials$ab) & !(x1 %in% antimicrobials$name)) + })]) > 0)) { + names(x)[which(str_detect(names(x), "[a-z]+") & !str_detect(names(x), "fullname") & !str_detect(names(x), "mo"))] <- str_to_upper(names(x)[which(str_detect(names(x), "[a-z]+") & !str_detect(names(x), "fullname") & !str_detect(names(x), "mo"))]) + names(x)[which(str_detect(names(x), "^ORG") & str_detect(names(x), "NAME$"))] <- "fullname" + names(x)[which(str_detect(names(x), "^ORG") & str_detect(names(x), "CODE$"))] <- "organism" + names(x)[which(str_detect(names(x), "^DATE") & str_detect(names(x), "CULTURE$") | str_detect(names(x), "^DATE") & str_detect(names(x), "SPECIMEN"))] <- "date" + names(x)[which(str_detect(names(x), "^EPISODE") | str_detect(names(x), "^PAT"))] <- "patient" + names(x)[which(str_detect(names(x), "^WARD"))] <- "ward" + names(x)[which(str_detect(names(x), "^HOSP"))] <- "hospital" + names(x)[which(str_detect(names(x), "DATE") & (str_detect(names(x), "REGISTRATION")))] <- "date" } - - #remove any duplicated variables - if (any(duplicated(names(x)))){ - x <- x[,-max(which(names(x) == names(x)[which(duplicated(names(x)))]))] + + # remove any duplicated variables + if (any(duplicated(names(x)))) { + x <- x[, -max(which(names(x) == names(x)[which(duplicated(names(x)))]))] } else { x <- x } - - x <- arrange(x,"date") #arrange in increasing date - - #if your dataset does not have an episode identifier - if (length(names(x)[sapply(names(x), function(x1){any(str_detect(x1,"patient"))})]) == 0){ - x$patient <- seq(1,nrow(x),1) + + x <- arrange(x, "date") # arrange in increasing date + + # if your dataset does not have an episode identifier + if (length(names(x)[sapply(names(x), function(x1) { + any(str_detect(x1, "patient")) + })]) == 0) { + x$patient <- seq(1, nrow(x), 1) } - - - - #find demographic data within dataset - micro_df_ind <- which(sapply(x,function(y) any(str_detect(y,"^S$"))|any(str_detect(y,"^R$"))|any(str_detect(y,"^SENSITIVE$"))| any(str_detect(y,"^RESISTANT$"))) == TRUE) - micro_df <- names(x)[micro_df_ind]#column names for susceptibility data - demograph_df <- names(x)[which(!names(x) %in% micro_df)] #column names for demographic or organism data - - #Dataset required to contain both mo and fullname - if ("mo" %in% names(x)){ + + + # find demographic data within dataset + micro_df_ind <- which(sapply(x, function(y) any(str_detect(y, "^S$")) | any(str_detect(y, "^R$")) | any(str_detect(y, "^SENSITIVE$")) | any(str_detect(y, "^RESISTANT$"))) == TRUE) + micro_df <- names(x)[micro_df_ind] # column names for susceptibility data + demograph_df <- names(x)[which(!names(x) %in% micro_df)] # column names for demographic or organism data + + # Dataset required to contain both mo and fullname + if ("mo" %in% names(x)) { # If mo included but not fullname: x <- x %>% mutate(mo = as.mo(mo), keep_synonyms = TRUE) test_isolates_pre <- x %>% - mutate(across(.cols = all_of(demograph_df),~ str_to_title(.x))) %>% - mutate(mo = toupper(mo)) - if ("fullname" %in% names(x)){ - test_isolates <- test_isolates_pre[,c(demograph_df)] + mutate(across(.cols = all_of(demograph_df), ~ str_to_title(.x))) %>% + mutate(mo = toupper(mo)) + if ("fullname" %in% names(x)) { + test_isolates <- test_isolates_pre[, c(demograph_df)] } else { - test_isolates <- left_join_microorganisms(x = test_isolates_pre, by = "mo")[,c(demograph_df,"fullname")] #add in fullname + test_isolates <- left_join_microorganisms(x = test_isolates_pre, by = "mo")[, c(demograph_df, "fullname")] # add in fullname } - test_isolates <- test_isolates %>% #remove bugs to be excluded and final style conversions + test_isolates <- test_isolates %>% # remove bugs to be excluded and final style conversions filter(!(fullname %in% exclude)) %>% mutate(fullname = str_to_title(fullname)) - - rx_abo <- x %>% #from orginal dataset - filter(mo %in% test_isolates$mo) %>% #select only those bugs that should be included - distinct(date, patient, mo, .keep_all = TRUE) #ensure that there are no duplications - rx_regs <- names(rx_abo)[!(names(rx_abo) %in% demograph_df)] #extract susceptibility data - - if (any(str_detect(rx_regs,fixed("+")))){#if combinations have already been accounted for in susceptibility profiling + + rx_abo <- x %>% # from orginal dataset + filter(mo %in% test_isolates$mo) %>% # select only those bugs that should be included + distinct(date, patient, mo, .keep_all = TRUE) # ensure that there are no duplications + rx_regs <- names(rx_abo)[!(names(rx_abo) %in% demograph_df)] # extract susceptibility data + + if (any(str_detect(rx_regs, fixed("+")))) { # if combinations have already been accounted for in susceptibility profiling rx_regs_ind <- which(!(names(rx_abo) %in% demograph_df)) # extract susceptibility data - rx_regs_comb <- rx_regs[str_detect(rx_regs,fixed("+"))]# find combination regimen variables + rx_regs_comb <- rx_regs[str_detect(rx_regs, fixed("+"))] # find combination regimen variables rx_regs_comb_x_total <- data.frame() - - for (jj in 1:length(rx_regs_comb)){ # for every combination regimen - rx_regs_comb_sep <- str_remove_all(unlist(str_split(rx_regs_comb[jj], fixed("+"))),"[[:punct:]]") #extract separate antimicrobials included in regimen + + for (jj in 1:length(rx_regs_comb)) { # for every combination regimen + rx_regs_comb_sep <- str_remove_all(unlist(str_split(rx_regs_comb[jj], fixed("+"))), "[[:punct:]]") # extract separate antimicrobials included in regimen rx_regs_comb_pre <- data.frame(x = rep(0, times = 2), y = rep(0, times = 2)) - for (ii in 1:length(rx_regs_comb_sep)){ #for every antibiotic in the combination regimen - rx_regs_comb_name <- ab_name(rx_regs_comb_sep[ii], only_first = TRUE) #extract their antibiotic name - rx_regs_comb_ab <- as.ab(rx_regs_comb_name) #convert antibiotic name to abbreviation - rx_regs_comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% rx_regs_comb_ab)[,c("ab","name")] - rx_regs_comb_pre[ii,] <- rx_regs_comb_antibiotic + for (ii in 1:length(rx_regs_comb_sep)) { # for every antibiotic in the combination regimen + rx_regs_comb_name <- ab_name(rx_regs_comb_sep[ii], only_first = TRUE) # extract their antibiotic name + rx_regs_comb_ab <- as.ab(rx_regs_comb_name) # convert antibiotic name to abbreviation + rx_regs_comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% rx_regs_comb_ab)[, c("ab", "name")] + rx_regs_comb_pre[ii, ] <- rx_regs_comb_antibiotic names(rx_regs_comb_pre) <- names(rx_regs_comb_antibiotic) } rx_regs_comb_x <- pivot_wider(rx_regs_comb_pre, names_from = name, values_from = ab) rx_regs_comb_x_set <- rx_regs_comb_x %>% - mutate(namerx = str_flatten(str_c(names(rx_regs_comb_x), collapse = "+"))) #generate new antimicrobials table - rx_regs_new_ab <- sapply(X = as.character(1:(ncol(rx_regs_comb_x_set)-1)), FUN = function(x){str_flatten(str_c("ab",x))}) - names(rx_regs_comb_x_set)[1:(ncol(rx_regs_comb_x_set)-1)] <- rx_regs_new_ab - names(rx_regs_comb_x_set)[ncol(rx_regs_comb_x_set)]<- "namerx" - rx_regs_comb_x_total <-rbind.data.frame(rx_regs_comb_x_total, rx_regs_comb_x_set, row.names = NULL) + mutate(namerx = str_flatten(str_c(names(rx_regs_comb_x), collapse = "+"))) # generate new antimicrobials table + rx_regs_new_ab <- sapply(X = as.character(1:(ncol(rx_regs_comb_x_set) - 1)), FUN = function(x) { + str_flatten(str_c("ab", x)) + }) + names(rx_regs_comb_x_set)[1:(ncol(rx_regs_comb_x_set) - 1)] <- rx_regs_new_ab + names(rx_regs_comb_x_set)[ncol(rx_regs_comb_x_set)] <- "namerx" + rx_regs_comb_x_total <- rbind.data.frame(rx_regs_comb_x_total, rx_regs_comb_x_set, row.names = NULL) } - - rx_regs_comb_rename <- select(rx_regs_comb_x_total,starts_with("ab")) #extract all antibiotic abbreviations - #regenerate variable names for dataset + + rx_regs_comb_rename <- select(rx_regs_comb_x_total, starts_with("ab")) # extract all antibiotic abbreviations + # regenerate variable names for dataset suppressWarnings({ - for (k in 1:nrow(rx_regs_comb_x_total)){ - rx_regs_comb_rename_2 <- unlist(rx_regs_comb_rename[k,]) + for (k in 1:nrow(rx_regs_comb_x_total)) { + rx_regs_comb_rename_2 <- unlist(rx_regs_comb_rename[k, ]) rx_regs_comb_rename_vec <- str_flatten(str_c(rx_regs_comb_rename_2, collapse = "+")) names(rx_abo)[names(rx_abo) == rx_regs_comb[k]] <- rx_regs_comb_rename_vec - } }) - names(rx_abo)[!(str_detect(names(rx_abo),fixed("+"))) & !(names(rx_abo) %in% demograph_df)] <- - as.ab(names(rx_abo)[!(str_detect(names(rx_abo),fixed("+"))) & !(names(rx_abo) %in% demograph_df)]) - names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))],".x") - rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "mo")[,c(names(rx_abo),"fullname")] - + names(rx_abo)[!(str_detect(names(rx_abo), fixed("+"))) & !(names(rx_abo) %in% demograph_df)] <- + as.ab(names(rx_abo)[!(str_detect(names(rx_abo), fixed("+"))) & !(names(rx_abo) %in% demograph_df)]) + names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))], ".x") + rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "mo")[, c(names(rx_abo), "fullname")] } else { names(rx_abo)[!(names(rx_abo) %in% demograph_df)] <- as.ab(names(rx_abo)[!(names(rx_abo) %in% demograph_df)]) - names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))],".x") - rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "mo")[,c(names(rx_abo),"fullname")] - - } + names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))], ".x") + rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "mo")[, c(names(rx_abo), "fullname")] + } } else { - #If fullname included in dataset: + # If fullname included in dataset: test_isolates_pre <- x %>% - mutate(across(.cols = all_of(demograph_df),~ str_to_title(.x))) %>% + mutate(across(.cols = all_of(demograph_df), ~ str_to_title(.x))) %>% filter(!(fullname %in% unique(exclude))) print(nrow(test_isolates_pre)) - test_isolates <- left_join_microorganisms(x = test_isolates_pre, by = "fullname")[,c(demograph_df,"mo")] - test_isolates$fullname <- str_to_title(test_isolates$fullname) + test_isolates <- left_join_microorganisms(x = test_isolates_pre, by = "fullname")[, c(demograph_df, "mo")] + test_isolates$fullname <- str_to_title(test_isolates$fullname) test_isolates$mo <- toupper(test_isolates$mo) - + rx_abo <- x %>% mutate(fullname = str_to_title(fullname)) %>% - filter(!(fullname %in% unique(exclude))) %>% + filter(!(fullname %in% unique(exclude))) %>% distinct(date, patient, fullname, .keep_all = TRUE) - rx_regs <- names(rx_abo)[!(names(rx_abo) %in% demograph_df)] #susceptibility data - - if (any(str_detect(rx_regs,fixed("+")))){#if combinations have already been accounted for + rx_regs <- names(rx_abo)[!(names(rx_abo) %in% demograph_df)] # susceptibility data + + if (any(str_detect(rx_regs, fixed("+")))) { # if combinations have already been accounted for rx_regs_ind <- which(!(names(rx_abo) %in% demograph_df)) - rx_regs_comb <- rx_regs[str_detect(rx_regs,fixed("+"))] + rx_regs_comb <- rx_regs[str_detect(rx_regs, fixed("+"))] rx_regs_comb_x_total <- data.frame() - - for (jj in 1:length(rx_regs_comb)){ - rx_regs_comb_sep <- str_remove_all(unlist(str_split(rx_regs_comb[jj], fixed("+"))),"[[:punct:]]") + + for (jj in 1:length(rx_regs_comb)) { + rx_regs_comb_sep <- str_remove_all(unlist(str_split(rx_regs_comb[jj], fixed("+"))), "[[:punct:]]") rx_regs_comb_pre <- data.frame(x = rep(0, times = 2), y = rep(0, times = 2)) - for (ii in 1:length(rx_regs_comb_sep)){ + for (ii in 1:length(rx_regs_comb_sep)) { rx_regs_comb_name <- ab_name(rx_regs_comb_sep[ii], only_first = TRUE) rx_regs_comb_ab <- as.ab(rx_regs_comb_name) - rx_regs_comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% rx_regs_comb_ab)[,c("ab","name")] - rx_regs_comb_pre[ii,] <- rx_regs_comb_antibiotic + rx_regs_comb_antibiotic <- subset(antimicrobials, antimicrobials$ab %in% rx_regs_comb_ab)[, c("ab", "name")] + rx_regs_comb_pre[ii, ] <- rx_regs_comb_antibiotic names(rx_regs_comb_pre) <- names(rx_regs_comb_antibiotic) } rx_regs_comb_x <- pivot_wider(rx_regs_comb_pre, names_from = name, values_from = ab) rx_regs_comb_x_set <- rx_regs_comb_x %>% mutate(namerx = str_flatten(str_c(names(rx_regs_comb_x), collapse = "+"))) - rx_regs_new_ab <- sapply(X = as.character(1:(ncol(rx_regs_comb_x_set)-1)), FUN = function(x){str_flatten(str_c("ab",x))}) - names(rx_regs_comb_x_set)[1:(ncol(rx_regs_comb_x_set)-1)] <- rx_regs_new_ab - names(rx_regs_comb_x_set)[ncol(rx_regs_comb_x_set)]<- "namerx" - rx_regs_comb_x_total <-rbind.data.frame(rx_regs_comb_x_total, rx_regs_comb_x_set, row.names = NULL) + rx_regs_new_ab <- sapply(X = as.character(1:(ncol(rx_regs_comb_x_set) - 1)), FUN = function(x) { + str_flatten(str_c("ab", x)) + }) + names(rx_regs_comb_x_set)[1:(ncol(rx_regs_comb_x_set) - 1)] <- rx_regs_new_ab + names(rx_regs_comb_x_set)[ncol(rx_regs_comb_x_set)] <- "namerx" + rx_regs_comb_x_total <- rbind.data.frame(rx_regs_comb_x_total, rx_regs_comb_x_set, row.names = NULL) } - - rx_regs_comb_rename <- select(rx_regs_comb_x_total,starts_with("ab")) - + + rx_regs_comb_rename <- select(rx_regs_comb_x_total, starts_with("ab")) + suppressWarnings({ - for (k in 1:nrow(rx_regs_comb_x_total)){ - rx_regs_comb_rename_2 <- unlist(rx_regs_comb_rename[k,]) + for (k in 1:nrow(rx_regs_comb_x_total)) { + rx_regs_comb_rename_2 <- unlist(rx_regs_comb_rename[k, ]) rx_regs_comb_rename_vec <- str_flatten(str_c(rx_regs_comb_rename_2, collapse = "+")) names(rx_abo)[names(rx_abo) == rx_regs_comb[k]] <- rx_regs_comb_rename_vec } }) - names(rx_abo)[!(str_detect(names(rx_abo),fixed("+"))) & !(names(rx_abo) %in% demograph_df)] <- - as.ab(names(rx_abo)[!(str_detect(names(rx_abo),fixed("+"))) & !(names(rx_abo) %in% demograph_df)]) - names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))],".x") - rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "fullname")[,c(names(rx_abo),"mo")] - + names(rx_abo)[!(str_detect(names(rx_abo), fixed("+"))) & !(names(rx_abo) %in% demograph_df)] <- + as.ab(names(rx_abo)[!(str_detect(names(rx_abo), fixed("+"))) & !(names(rx_abo) %in% demograph_df)]) + names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))], ".x") + rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "fullname")[, c(names(rx_abo), "mo")] } else { names(rx_abo)[!(names(rx_abo) %in% demograph_df)] <- as.ab(names(rx_abo)[!(names(rx_abo) %in% demograph_df)]) - names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))],".x") - rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "fullname")[,c(names(rx_abo),"mo")] - - } + names(rx_abo)[which(duplicated(names(rx_abo)))] <- str_c(names(rx_abo)[which(duplicated(names(rx_abo)))], ".x") + rx_abo_upd <- left_join_microorganisms(x = rx_abo, by = "fullname")[, c(names(rx_abo), "mo")] + } } - - + + rx_abo_upd$fullname <- str_to_title(rx_abo_upd$fullname) - - - #generate new dataset with demographic information and only antimicrobials under investigation - rx_1 <- select(antibiotic_rx,starts_with("ab")) #all antimicrobials + + + # generate new dataset with demographic information and only antimicrobials under investigation + rx_1 <- select(antibiotic_rx, starts_with("ab")) # all antimicrobials rx_df <- test_isolates - var_date <- names(rx_df)[str_detect(names(rx_df),"date")] + var_date <- names(rx_df)[str_detect(names(rx_df), "date")] rx_df[[var_date]] <- ymd(rx_df[[var_date]]) - var_patient <- names(rx_df)[str_detect(names(rx_df),"patient")] + var_patient <- names(rx_df)[str_detect(names(rx_df), "patient")] rx_df[[var_patient]] <- toupper(rx_df[[var_patient]]) - + names_vec <- vector() suppressWarnings({ - for (k in 1:nrow(rx_1)){ - rx_2 <- unlist(rx_1[k,][!is.na(str_extract(rx_1[k,],"[A-Z]+"))]) - if ((length(rx_2) > 1) & !(str_flatten(str_c(rx_2, collapse = "+")) %in% names(rx_abo_upd))){ # if more than one antibiotic in regimen and susceptibility of the combination not considered - rx_df <- merge(x = rx_df, y = rx_abo_upd[,c("date","patient","mo",rx_2)], by = c("date","patient","mo")) #then bind separate antimicrobials in regimen - #if prior individual antibiotic was added, suffix will be generated - for (jj in 1:length(rx_2)){ - if (any(str_detect(names(rx_df),paste0(rx_2[jj],".y")))){ - rx_df[[paste0(rx_2[jj],".y")]] <- str_trim(rx_df[[paste0(rx_2[jj],".y")]]) + for (k in 1:nrow(rx_1)) { + rx_2 <- unlist(rx_1[k, ][!is.na(str_extract(rx_1[k, ], "[A-Z]+"))]) + if ((length(rx_2) > 1) & !(str_flatten(str_c(rx_2, collapse = "+")) %in% names(rx_abo_upd))) { # if more than one antibiotic in regimen and susceptibility of the combination not considered + rx_df <- merge(x = rx_df, y = rx_abo_upd[, c("date", "patient", "mo", rx_2)], by = c("date", "patient", "mo")) # then bind separate antimicrobials in regimen + # if prior individual antibiotic was added, suffix will be generated + for (jj in 1:length(rx_2)) { + if (any(str_detect(names(rx_df), paste0(rx_2[jj], ".y")))) { + rx_df[[paste0(rx_2[jj], ".y")]] <- str_trim(rx_df[[paste0(rx_2[jj], ".y")]]) } else { rx_df[[rx_2[jj]]] <- str_trim(rx_df[[rx_2[jj]]]) } } - names(rx_df)[names(rx_df) %in% rx_2] <- str_flatten(str_c(rx_2, collapse = "+")) #rename both variables as the combination - names(rx_df)[names(rx_df) %in% paste0(rx_2,".y")] <- str_flatten(str_c(rx_2, collapse = "+")) - names(rx_df)[names(rx_df) %in% paste0(rx_2,".x")] <- rx_2 + names(rx_df)[names(rx_df) %in% rx_2] <- str_flatten(str_c(rx_2, collapse = "+")) # rename both variables as the combination + names(rx_df)[names(rx_df) %in% paste0(rx_2, ".y")] <- str_flatten(str_c(rx_2, collapse = "+")) + names(rx_df)[names(rx_df) %in% paste0(rx_2, ".x")] <- rx_2 names_vec[k] <- str_flatten(str_c(rx_2, collapse = "+")) - } else if ((length(rx_2) > 1) & (str_flatten(str_c(rx_2, collapse = "+")) %in% names(rx_abo_upd))){#if more than one antibiotic in regimen and susceptibility of the combination considered - rx_df <- merge(x = rx_df, y = rx_abo_upd[,c("date","patient","mo",str_flatten(str_c(rx_2, collapse = "+")))], by = c("date","patient","mo")) + } else if ((length(rx_2) > 1) & (str_flatten(str_c(rx_2, collapse = "+")) %in% names(rx_abo_upd))) { # if more than one antibiotic in regimen and susceptibility of the combination considered + rx_df <- merge(x = rx_df, y = rx_abo_upd[, c("date", "patient", "mo", str_flatten(str_c(rx_2, collapse = "+")))], by = c("date", "patient", "mo")) rx_df[[str_flatten(str_c(rx_2, collapse = "+"))]] <- str_trim(rx_df[[str_flatten(str_c(rx_2, collapse = "+"))]]) names_vec[k] <- str_flatten(str_c(rx_2, collapse = "+")) } else { - rx_df <- merge(x = rx_df, y = rx_abo_upd[,c("date","patient","mo",rx_2)], by = c("date","patient","mo")) #then bind separate antimicrobials in regimen + rx_df <- merge(x = rx_df, y = rx_abo_upd[, c("date", "patient", "mo", rx_2)], by = c("date", "patient", "mo")) # then bind separate antimicrobials in regimen rx_df[[rx_2]] <- str_trim(rx_df[[rx_2]]) - names(rx_df)[names(rx_df) %in% rx_2] <- str_flatten(str_c(rx_2, collapse = "+")) #rename both variables as the combination - names(rx_df)[names(rx_df) %in% paste0(rx_2,".y")] <- str_flatten(str_c(rx_2, collapse = "+")) - names(rx_df)[names(rx_df) %in% paste0(rx_2,".x")] <- rx_2 + names(rx_df)[names(rx_df) %in% rx_2] <- str_flatten(str_c(rx_2, collapse = "+")) # rename both variables as the combination + names(rx_df)[names(rx_df) %in% paste0(rx_2, ".y")] <- str_flatten(str_c(rx_2, collapse = "+")) + names(rx_df)[names(rx_df) %in% paste0(rx_2, ".x")] <- rx_2 names_vec[k] <- rx_2 } } }) - - if (any(is.na(rx_df$mo))){ + + if (any(is.na(rx_df$mo))) { missing1 <- which(is.na(rx_df$mo)) - if (rx_df[missing1,"fullname"] != "NA"){ - rx_df[missing1,"mo"] <- as.mo(rx_df[missing1,]$fullname, keep_synonyms = TRUE) + if (rx_df[missing1, "fullname"] != "NA") { + rx_df[missing1, "mo"] <- as.mo(rx_df[missing1, ]$fullname, keep_synonyms = TRUE) } else { - rx_df <- rx_df[-missing1,] + rx_df <- rx_df[-missing1, ] } } - - - #Generate adapted first isolate dataset on antimicrobial level - if (isolate_first == "yes"){ #if dataset has already been deduplicated for first isolates + + + # Generate adapted first isolate dataset on antimicrobial level + if (isolate_first == "yes") { # if dataset has already been deduplicated for first isolates test_isolates_first <- rx_df } else { - test_isolates_first <- rx_df[first_isolate(x = rx_df, - method = "e", - episode_days = 14, - col_date = var_date, - col_mo = "mo") == TRUE,] #Only include isolates for microorganisms < 14 days apart + test_isolates_first <- rx_df[first_isolate( + x = rx_df, + method = "e", + episode_days = 14, + col_date = var_date, + col_mo = "mo" + ) == TRUE, ] # Only include isolates for microorganisms < 14 days apart } - - - + + test_isolates_first_upd <<- test_isolates_first %>% - pivot_longer(cols = names(test_isolates_first)[which(sapply(test_isolates_first, - function(y) any(str_detect(y,"^S$"))|any(str_detect(y,"^R$"))|any(str_detect(y,"^SENSITIVE$"))| any(str_detect(y,"^RESISTANT$"))))], - names_to = "keyantimicrobials", values_to = "antibiogram") %>% #combine all regimens and antibiograms into single variables - mutate(antibiogram = ifelse(antibiogram == "I"|antibiogram == "INTERMEDIATE",susceptible_I,antibiogram)) %>% #simplify susceptibility outputs - 90-60 rule - mutate(antibiogram = ifelse(antibiogram == "SENSITIVE","S", antibiogram)) %>% - mutate(antibiogram = ifelse(antibiogram == "RESISTANT","R",antibiogram)) %>% - mutate(antibiogram = factor(antibiogram), - keyantimicrobials = as.factor(keyantimicrobials)) - - - #recode susceptibility of regimes if not already in the dataset: if any S, then all S; if any NA and not S, then combination NA; if all R, then combination R. + pivot_longer( + cols = names(test_isolates_first)[which(sapply( + test_isolates_first, + function(y) any(str_detect(y, "^S$")) | any(str_detect(y, "^R$")) | any(str_detect(y, "^SENSITIVE$")) | any(str_detect(y, "^RESISTANT$")) + ))], + names_to = "keyantimicrobials", values_to = "antibiogram" + ) %>% # combine all regimens and antibiograms into single variables + mutate(antibiogram = ifelse(antibiogram == "I" | antibiogram == "INTERMEDIATE", susceptible_I, antibiogram)) %>% # simplify susceptibility outputs - 90-60 rule + mutate(antibiogram = ifelse(antibiogram == "SENSITIVE", "S", antibiogram)) %>% + mutate(antibiogram = ifelse(antibiogram == "RESISTANT", "R", antibiogram)) %>% + mutate( + antibiogram = factor(antibiogram), + keyantimicrobials = as.factor(keyantimicrobials) + ) + + + # recode susceptibility of regimes if not already in the dataset: if any S, then all S; if any NA and not S, then combination NA; if all R, then combination R. print("-------------------") print(test_isolates_first_upd) print(names(test_isolates_first_upd)) print(duplicated(names(test_isolates_first_upd))) print(any(duplicated(names(test_isolates_first_upd)))) print("-------------------") - - if (any(duplicated(names(test_isolates_first_upd)))){ + + if (any(duplicated(names(test_isolates_first_upd)))) { test_isolates_temp3 <- test_isolates_upd %>% mutate(antibiogram = as.character(antibiogram)) %>% - mutate(antibiogram = replace_na(antibiogram,"U")) %>% + mutate(antibiogram = replace_na(antibiogram, "U")) %>% mutate(combiantibiogram = "A") test_isolates_temp4 <- test_isolates_temp3 %>% - group_by(date,patient,mo,keyantimicrobials) %>% - mutate(combiantibiogram = ifelse(any(str_detect(antibiogram,"S")), "S",combiantibiogram)) %>% - mutate(combiantibiogram = ifelse(any(str_detect(antibiogram,"U")) & !any(str_detect(antibiogram,"S")), "U",combiantibiogram)) %>% - mutate(combiantibiogram = ifelse(!any(str_detect(antibiogram,"U")) & !any(str_detect(antibiogram,"S")), "R", combiantibiogram)) %>% - distinct(date,patient,mo,.keep_all = TRUE) %>% + group_by(date, patient, mo, keyantimicrobials) %>% + mutate(combiantibiogram = ifelse(any(str_detect(antibiogram, "S")), "S", combiantibiogram)) %>% + mutate(combiantibiogram = ifelse(any(str_detect(antibiogram, "U")) & !any(str_detect(antibiogram, "S")), "U", combiantibiogram)) %>% + mutate(combiantibiogram = ifelse(!any(str_detect(antibiogram, "U")) & !any(str_detect(antibiogram, "S")), "R", combiantibiogram)) %>% + distinct(date, patient, mo, .keep_all = TRUE) %>% ungroup() %>% select(-antibiogram) %>% rename("antibiogram" = "combiantibiogram") %>% - mutate(antibiogram = ifelse(antibiogram == "U",NA_character_, antibiogram)) %>% + mutate(antibiogram = ifelse(antibiogram == "U", NA_character_, antibiogram)) %>% mutate(antibiogram = as.factor(antibiogram)) - - test_isolates_temp5 <- test_isolates_upd %>% #isolate the single regimens - group_by(date,patient,mo) %>% - filter(!str_detect(keyantimicrobials,fixed("+"))) %>% + + test_isolates_temp5 <- test_isolates_upd %>% # isolate the single regimens + group_by(date, patient, mo) %>% + filter(!str_detect(keyantimicrobials, fixed("+"))) %>% ungroup() - - test_isolates_final <- bind_rows(test_isolates_temp4,test_isolates_temp5) + + test_isolates_final <- bind_rows(test_isolates_temp4, test_isolates_temp5) test_isolates_final <- test_isolates_final %>% - arrange(patient) - + arrange(patient) } else { test_isolates_final <- test_isolates_first_upd - } - + } + test_isolates_final <<- test_isolates_final - + test_isolates_final$fullname <- as.factor(test_isolates_final$fullname) test_isolates_final$date <- ymd(test_isolates_final$date) - - - #Calculate pathogen incidence and susceptible by regimen - mo_name_total <- test_isolates_final %>% #deduplicate pathogen list - filter(!duplicated(mo)) %>% - select(mo,fullname) - - organism_count_total <- test_isolates_final %>% #count the number of occurrences of pathogens within first isolates set + + + # Calculate pathogen incidence and susceptible by regimen + mo_name_total <- test_isolates_final %>% # deduplicate pathogen list + filter(!duplicated(mo)) %>% + select(mo, fullname) + + organism_count_total <- test_isolates_final %>% # count the number of occurrences of pathogens within first isolates set group_by(mo) %>% - distinct(date,patient,.keep_all = TRUE) %>% + distinct(date, patient, .keep_all = TRUE) %>% count(mo) %>% ungroup() %>% arrange(desc(n)) - + organism_count_total_upd <- organism_count_total %>% mutate(mo = as.character(mo)) %>% - mutate(mo = ifelse(str_detect(mo, "^B_KLBSL"),"B_KLBSL",mo)) %>% - mutate(mo = as.factor(mo)) + mutate(mo = ifelse(str_detect(mo, "^B_KLBSL"), "B_KLBSL", mo)) %>% + mutate(mo = as.factor(mo)) klebsiella_combine <- organism_count_total_upd %>% filter(mo == "B_KLBSL") %>% summarise(n_kleb = sum(n)) organism_count_total_upd <- organism_count_total_upd %>% mutate(n = ifelse(mo == "B_KLBSL", klebsiella_combine$n_kleb, n)) %>% - distinct(mo,.keep_all = TRUE) + distinct(mo, .keep_all = TRUE) organism_count_total_upd <- rename(organism_count_total_upd, - "n_inc" = "n") - - regimen_organism <- slice_max(organism_count_total_upd, order_by = n_inc, n = pathogen_in) #choose only the top n number of occurring pathogens to consider in WISCA - total_pathogens <- sum(regimen_organism$n_inc) #total number of occurrences of pathogens - print(paste0("Total pathogens = ",sum(organism_count_total_upd$n_inc))) + "n_inc" = "n" + ) + + regimen_organism <- slice_max(organism_count_total_upd, order_by = n_inc, n = pathogen_in) # choose only the top n number of occurring pathogens to consider in WISCA + total_pathogens <- sum(regimen_organism$n_inc) # total number of occurrences of pathogens + print(paste0("Total pathogens = ", sum(organism_count_total_upd$n_inc))) print(paste0("Number of pathogens contributing to WISCA = ", total_pathogens)) - print(paste0("Percentage of pathogens contributing to WISCA = ",round((total_pathogens/sum(organism_count_total_upd$n_inc))*100,0),"%")) - + print(paste0("Percentage of pathogens contributing to WISCA = ", round((total_pathogens / sum(organism_count_total_upd$n_inc)) * 100, 0), "%")) + test_isolates_final <- test_isolates_final %>% mutate(mo = as.character(mo)) %>% - mutate(mo = ifelse(str_detect(mo, "^B_KLBSL"), "B_KLBSL",mo)) %>% #combine all Klebsiella species together + mutate(mo = ifelse(str_detect(mo, "^B_KLBSL"), "B_KLBSL", mo)) %>% # combine all Klebsiella species together mutate(mo = as.factor(mo)) %>% mutate(fullname = as.character(fullname)) %>% - mutate(fullname = ifelse(str_detect(fullname, "Klebsiella"), "Klebsiella sp.", fullname)) %>% #rename all Klebsiella species + mutate(fullname = ifelse(str_detect(fullname, "Klebsiella"), "Klebsiella sp.", fullname)) %>% # rename all Klebsiella species mutate(fullname = as.factor(fullname)) - - - #for each antimicrobial regimen calculate the pathogen incidence and sensitivity + + + # for each antimicrobial regimen calculate the pathogen incidence and sensitivity params_set_out <- data.frame() - - if (analysis %in% names(test_isolates_final)){ - - #if there is a facility/ward/category level analysis required - for (j in 1:length(unique(test_isolates_final[[analysis]]))){#for every facility/ward/category - infection_episodes <- length(unique(test_isolates_final[test_isolates_final[[analysis]] == unique(test_isolates_final[[analysis]])[j],]$patient)) - - if (infection_episodes < 100){#must have at least 100 infection episodes for WISCA - warning(paste0("Number of infection episodes = ",infection_episodes,". Minimum sample size for WISCA is 100 infection episodes. Below minimum sample size will produce inaccurate coverage estimates.")) - + + if (analysis %in% names(test_isolates_final)) { + # if there is a facility/ward/category level analysis required + for (j in 1:length(unique(test_isolates_final[[analysis]]))) { # for every facility/ward/category + infection_episodes <- length(unique(test_isolates_final[test_isolates_final[[analysis]] == unique(test_isolates_final[[analysis]])[j], ]$patient)) + + if (infection_episodes < 100) { # must have at least 100 infection episodes for WISCA + warning(paste0("Number of infection episodes = ", infection_episodes, ". Minimum sample size for WISCA is 100 infection episodes. Below minimum sample size will produce inaccurate coverage estimates.")) } else { - print(paste0("Number of infection episodes = ",infection_episodes," > 100 minimum required. Sample size adequate for WISCA analysis.")) + print(paste0("Number of infection episodes = ", infection_episodes, " > 100 minimum required. Sample size adequate for WISCA analysis.")) } - - + + y_site <- subset(test_isolates_final, test_isolates_final[[analysis]] == unique(test_isolates_final[[analysis]])[j]) - - mo_name_sites <- y_site %>% #deduplicate - filter(!duplicated(mo)) %>% - select(mo,fullname,all_of(analysis)) - - organism_count_regimen <- y_site %>% #count the number of organisms per regimen + + mo_name_sites <- y_site %>% # deduplicate + filter(!duplicated(mo)) %>% + select(mo, fullname, all_of(analysis)) + + organism_count_regimen <- y_site %>% # count the number of organisms per regimen filter(mo %in% regimen_organism$mo) %>% group_by(mo) %>% distinct(date, patient, .keep_all = TRUE) %>% count(mo) %>% ungroup() %>% arrange(desc(n)) - - - organism_count_df <- left_join(x = organism_count_regimen, y = mo_name_sites, by = "mo")[,c("mo",analysis,"fullname","n")] - + + + organism_count_df <- left_join(x = organism_count_regimen, y = mo_name_sites, by = "mo")[, c("mo", analysis, "fullname", "n")] + total_sites_pathogens <- sum(organism_count_df$n) - print(paste0("Number of pathogens contributing to WISCA in ",unique(test_isolates_final[[analysis]])[j]," = ", total_sites_pathogens)) - - organism_count_df <- organism_count_df %>% - mutate(prop.n = n/total_sites_pathogens) #pathogen incidence = number of pathogen occurring/total pathogen occurrences included in WISCA - # - - susceptible_regimen_tested <<- y_site %>% #determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested + print(paste0("Number of pathogens contributing to WISCA in ", unique(test_isolates_final[[analysis]])[j], " = ", total_sites_pathogens)) + + organism_count_df <- organism_count_df %>% + mutate(prop.n = n / total_sites_pathogens) # pathogen incidence = number of pathogen occurring/total pathogen occurrences included in WISCA + # + + susceptible_regimen_tested <<- y_site %>% # determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested filter(mo %in% regimen_organism$mo) %>% filter(!is.na(antibiogram)) %>% - group_by(mo,keyantimicrobials) %>% + group_by(mo, keyantimicrobials) %>% distinct(date, patient, .keep_all = TRUE) %>% count(antibiogram) %>% ungroup() %>% - arrange(desc(n)) - susceptible_regimen_ntested <<- y_site %>% #determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested + arrange(desc(n)) + susceptible_regimen_ntested <<- y_site %>% # determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested filter(mo %in% regimen_organism$mo) %>% - filter(is.na(antibiogram)) %>% #Not tested - group_by(mo,keyantimicrobials) %>% + filter(is.na(antibiogram)) %>% # Not tested + group_by(mo, keyantimicrobials) %>% distinct(date, patient, .keep_all = TRUE) %>% count(antibiogram) %>% ungroup() %>% - arrange(desc(n)) - - organism_count_df_tested <<- susceptible_regimen_tested %>% - group_by(mo,keyantimicrobials) %>% + arrange(desc(n)) + + organism_count_df_tested <<- susceptible_regimen_tested %>% + group_by(mo, keyantimicrobials) %>% summarise(tested_n = sum(n), .groups = "keep") %>% ungroup() - organism_count_df_ntested <<- susceptible_regimen_ntested %>% - group_by(mo,keyantimicrobials) %>% + organism_count_df_ntested <<- susceptible_regimen_ntested %>% + group_by(mo, keyantimicrobials) %>% summarise(ntested_n = sum(n), .groups = "keep") %>% ungroup() - - organism_count_df_combine_tested <- full_join(x = organism_count_df_tested, y = organism_count_df_ntested, by = c("mo","keyantimicrobials")) #combine tested and nontested + + organism_count_df_combine_tested <- full_join(x = organism_count_df_tested, y = organism_count_df_ntested, by = c("mo", "keyantimicrobials")) # combine tested and nontested organism_count_df_combine_tested <- organism_count_df_combine_tested %>% - mutate(tested_n = ifelse(is.na(tested_n), 0 ,tested_n), - ntested_n = ifelse(is.na(ntested_n), 0, ntested_n)) - - organism_count_df_regimen_tested <- left_join(x = organism_count_df, y = organism_count_df_combine_tested[,c("mo","keyantimicrobials","tested_n")], by = "mo") - - #minimum number of tested isolates per ABO should be 30 for accurate coverage estimates + mutate( + tested_n = ifelse(is.na(tested_n), 0, tested_n), + ntested_n = ifelse(is.na(ntested_n), 0, ntested_n) + ) + + organism_count_df_regimen_tested <- left_join(x = organism_count_df, y = organism_count_df_combine_tested[, c("mo", "keyantimicrobials", "tested_n")], by = "mo") + + # minimum number of tested isolates per ABO should be 30 for accurate coverage estimates condition_tested <- organism_count_df_regimen_tested %>% - filter(tested_n <30) - if (nrow(condition_tested) > 0){ + filter(tested_n < 30) + if (nrow(condition_tested) > 0) { warning("Number of tested isolates per regimen should exceed 30. Coverage estimates will be inaccurate for following regimens") - View(condition_tested[c("mo","n","keyantimicrobials","tested_n")]) + View(condition_tested[c("mo", "n", "keyantimicrobials", "tested_n")]) } else { print("Number of tested isolates per regimen exceed minimum threshold of 30.") } - + organism_tested_sensitive <- susceptible_regimen_tested %>% filter(antibiogram == "S") - - if (nrow(susceptible_regimen_tested)> 0){ - sensitivity <- left_join(x = organism_count_df_regimen_tested, y = organism_tested_sensitive, by = c("mo","keyantimicrobials")) + + if (nrow(susceptible_regimen_tested) > 0) { + sensitivity <- left_join(x = organism_count_df_regimen_tested, y = organism_tested_sensitive, by = c("mo", "keyantimicrobials")) sensitivity <- sensitivity %>% rename("S_n" = "n.y") %>% - mutate(S_n = ifelse(is.na(S_n),0,S_n)) %>% + mutate(S_n = ifelse(is.na(S_n), 0, S_n)) %>% select(-antibiogram) - sensitivity <- sensitivity %>% #calculate proportion sensitive out of those tested and having results - mutate(prop.S = ifelse(tested_n != 0, S_n/tested_n, 0)) + sensitivity <- sensitivity %>% # calculate proportion sensitive out of those tested and having results + mutate(prop.S = ifelse(tested_n != 0, S_n / tested_n, 0)) sensitivity$perc.S <- sensitivity$prop.S * 100 - - params_set_out <- rbind.data.frame(params_set_out, sensitivity, row.names = NULL) - } else { + + params_set_out <- rbind.data.frame(params_set_out, sensitivity, row.names = NULL) + } else { next } } - + params_set_out$fullname <- droplevels(params_set_out$fullname) params_set_out[[names(test_isolates_final)[names(test_isolates_final) == analysis]]] <- as.factor(params_set_out[[names(test_isolates_final)[names(test_isolates_final) == analysis]]]) - - } else{ - + } else { infection_episodes <- length(unique(test_isolates_final$patient)) - if (infection_episodes < 100){#must have at least 100 infection episodes for WISCA - warning(paste0("Number of infection episodes = ",infection_episodes,". Minimum sample size for WISCA is 100 infection episodes. Below minimum sample size will produce inaccurate coverage estimates.")) - + if (infection_episodes < 100) { # must have at least 100 infection episodes for WISCA + warning(paste0("Number of infection episodes = ", infection_episodes, ". Minimum sample size for WISCA is 100 infection episodes. Below minimum sample size will produce inaccurate coverage estimates.")) } else { - print(paste0("Number of infection episodes = ",infection_episodes," > 100 minimum required. Sample size adequate for WISCA analysis.")) + print(paste0("Number of infection episodes = ", infection_episodes, " > 100 minimum required. Sample size adequate for WISCA analysis.")) } - - mo_name_sites <- test_isolates_final %>% #deduplicate - filter(!duplicated(mo)) %>% - select(mo,fullname) - - organism_count_path <- test_isolates_final %>% #count the number of organisms + + mo_name_sites <- test_isolates_final %>% # deduplicate + filter(!duplicated(mo)) %>% + select(mo, fullname) + + organism_count_path <- test_isolates_final %>% # count the number of organisms filter(mo %in% regimen_organism$mo) %>% group_by(mo) %>% distinct(date, patient, .keep_all = TRUE) %>% count(mo) %>% ungroup() %>% - arrange(desc(n)) - - organism_count_df <- left_join(x = organism_count_path, y = mo_name_sites, by = "mo")[,c("mo","fullname","n")] - + arrange(desc(n)) + + organism_count_df <- left_join(x = organism_count_path, y = mo_name_sites, by = "mo")[, c("mo", "fullname", "n")] + total_sites_pathogens <- sum(organism_count_df$n) - - - - organism_count_df <- organism_count_df %>% - mutate(prop.n = n/total_pathogens) #pathogen incidence = number of pathogen occurring/total pathogen occurrences included in WISCA - - - susceptible_regimen_tested <<- test_isolates_final %>% #determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested + + + organism_count_df <- organism_count_df %>% + mutate(prop.n = n / total_pathogens) # pathogen incidence = number of pathogen occurring/total pathogen occurrences included in WISCA + + + susceptible_regimen_tested <<- test_isolates_final %>% # determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested filter(mo %in% regimen_organism$mo) %>% - filter(!is.na(antibiogram)) %>% #Tested - group_by(mo,keyantimicrobials) %>% + filter(!is.na(antibiogram)) %>% # Tested + group_by(mo, keyantimicrobials) %>% distinct(date, patient, .keep_all = TRUE) %>% count(antibiogram) %>% ungroup() %>% - arrange(desc(n)) - susceptible_regimen_ntested <<- test_isolates_final %>% #determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested + arrange(desc(n)) + susceptible_regimen_ntested <<- test_isolates_final %>% # determine sensitivity profile - assume those with NA not tested. Calculate susceptibility only on tested filter(mo %in% regimen_organism$mo) %>% - filter(is.na(antibiogram)) %>% #Not tested - group_by(mo,keyantimicrobials) %>% + filter(is.na(antibiogram)) %>% # Not tested + group_by(mo, keyantimicrobials) %>% distinct(date, patient, .keep_all = TRUE) %>% count(antibiogram) %>% ungroup() %>% - arrange(desc(n)) - - organism_count_df_tested <<- susceptible_regimen_tested %>% - group_by(mo,keyantimicrobials) %>% + arrange(desc(n)) + + organism_count_df_tested <<- susceptible_regimen_tested %>% + group_by(mo, keyantimicrobials) %>% summarise(tested_n = sum(n), .groups = "keep") %>% ungroup() - organism_count_df_ntested <<- susceptible_regimen_ntested %>% - group_by(mo,keyantimicrobials) %>% + organism_count_df_ntested <<- susceptible_regimen_ntested %>% + group_by(mo, keyantimicrobials) %>% summarise(ntested_n = sum(n), .groups = "keep") %>% ungroup() - - organism_count_df_combine_tested <- full_join(x = organism_count_df_tested, y = organism_count_df_ntested, by = c("mo","keyantimicrobials")) #combine tested and nontested + + organism_count_df_combine_tested <- full_join(x = organism_count_df_tested, y = organism_count_df_ntested, by = c("mo", "keyantimicrobials")) # combine tested and nontested organism_count_df_combine_tested <- organism_count_df_combine_tested %>% - mutate(tested_n = ifelse(is.na(tested_n), 0 ,tested_n), - ntested_n = ifelse(is.na(ntested_n), 0, ntested_n)) - - organism_count_df_regimen_tested <- left_join(x = organism_count_df, y = organism_count_df_combine_tested[,c("mo","keyantimicrobials","tested_n")], by = "mo") - - #minimum number of tested isolates per ABO should be 30 for accurate coverage estimates + mutate( + tested_n = ifelse(is.na(tested_n), 0, tested_n), + ntested_n = ifelse(is.na(ntested_n), 0, ntested_n) + ) + + organism_count_df_regimen_tested <- left_join(x = organism_count_df, y = organism_count_df_combine_tested[, c("mo", "keyantimicrobials", "tested_n")], by = "mo") + + # minimum number of tested isolates per ABO should be 30 for accurate coverage estimates condition_tested <- organism_count_df_regimen_tested %>% - filter(tested_n <30) - if (nrow(condition_tested) > 0){ + filter(tested_n < 30) + if (nrow(condition_tested) > 0) { warning("Number of tested isolates per regimen should exceed 30. Coverage estimates will be inaccurate for these regimens.") - View(condition_tested[,c("mo","n","keyantimicrobials","tested_n")]) + View(condition_tested[, c("mo", "n", "keyantimicrobials", "tested_n")]) } else { print("Number of tested isolates per regimen exceed minimum threshold of 30.") } - - + + organism_tested_sensitive <- susceptible_regimen_tested %>% filter(antibiogram == "S") - - if (nrow(susceptible_regimen_tested)> 0){ - sensitivity <- left_join(x = organism_count_df_regimen_tested, y = organism_tested_sensitive, by = c("mo","keyantimicrobials")) + + if (nrow(susceptible_regimen_tested) > 0) { + sensitivity <- left_join(x = organism_count_df_regimen_tested, y = organism_tested_sensitive, by = c("mo", "keyantimicrobials")) sensitivity <- sensitivity %>% rename("S_n" = "n.y") %>% - mutate(S_n = ifelse(is.na(S_n),0,S_n)) %>% + mutate(S_n = ifelse(is.na(S_n), 0, S_n)) %>% select(-antibiogram) - sensitivity <- sensitivity %>% #calculate proportion sensitive out of those tested and having results - mutate(prop.S = ifelse(tested_n != 0, S_n/tested_n, 0)) + sensitivity <- sensitivity %>% # calculate proportion sensitive out of those tested and having results + mutate(prop.S = ifelse(tested_n != 0, S_n / tested_n, 0)) sensitivity$perc.S <- sensitivity$prop.S * 100 - - params_set_out <- rbind.data.frame(params_set_out, sensitivity, row.names = NULL) - + + params_set_out <- rbind.data.frame(params_set_out, sensitivity, row.names = NULL) } else { next } params_set_out$fullname <- droplevels(params_set_out$fullname) - } - + # if (infection_in != ""){ # params_set_out <- params_set_out %>% # mutate(infection_type = infection_in) %>% # mutate(infection_type_full = infection_full) %>% # rowwise() %>% # mutate(regimens = ifelse(str_detect(keyantimicrobials, fixed("+")), - # str_c(sapply(unlist(unique(str_split(keyantimicrobials, fixed("+")))),ab_name,USE.NAMES = FALSE), collapse = " + "), + # str_c(sapply(unlist(unique(str_split(keyantimicrobials, fixed("+")))),ab_name,USE.NAMES = FALSE), collapse = " + "), # ab_name(keyantimicrobials))) %>% # select(mo:prop.n, regimens,keyantimicrobials,tested_n:infection_type, infection_type_full) - # - # + # + # # } else { # print("here??") # params_set_out <- params_set_out %>% # mutate(infection_type_full = infection_full) %>% # rowwise() %>% # mutate(regimens = ifelse(str_detect(keyantimicrobials, fixed("+")), - # str_c(sapply(unlist(unique(str_split(keyantimicrobials, fixed("+")))),ab_name,USE.NAMES = FALSE), collapse = " + "), + # str_c(sapply(unlist(unique(str_split(keyantimicrobials, fixed("+")))),ab_name,USE.NAMES = FALSE), collapse = " + "), # ab_name(keyantimicrobials))) %>% # select(mo:prop.n, regimens,keyantimicrobials,tested_n:infection_type, infection_type_full) # } - - + return(params_set_out) - - } diff --git a/data-raw/wisca_test/wisca_AMR.R b/data-raw/wisca_test/wisca_AMR.R index d6f16c4eb..3db376b76 100644 --- a/data-raw/wisca_test/wisca_AMR.R +++ b/data-raw/wisca_test/wisca_AMR.R @@ -1,122 +1,128 @@ - -# Copyright (c) [2022] [Larisse Bolton (author), Aislinn Cook (contributor)] -# Adapted from: Bielicki JA, Sharland M, Heath PT, et al. Evaluation of the coverage of 3 antibiotic regimens for -# neonatal sepsis in the hospital setting across Asian countries. JAMA Netw Open. -# 2020:3(2):e1921124. doi:10.1001.jamanetworkopen.2019.21124 - - -wisca_funct <- function(wisca_in,analysis){ - -#function used to gerenrate coverage estimate parameters from input data -priors <- function(wisca_prior_in){ - #INCIDENCE - #prior Dirichlet (Gamma) parameters - A <- rep(1, times = length(unique(wisca_prior_in$mo))) - Multiplier <- rep(1, times = length(unique(wisca_prior_in$mo))) - multinomial_obs <- round((wisca_prior_in$n.x)*Multiplier,0) - - #posterior Dirichlet (Gamma) parameters - gamma_A <- A + multinomial_obs - B <- rep(1, times = length(gamma_A)) - - ##SUSCEPTIBILITY - #prior Beta parameters - beta_A <- rep(1, times = length(unique(wisca_prior_in$mo))) - beta_B <- rep(1, times = length(unique(wisca_prior_in$mo))) - - #Binomial distribution parameters - r <- wisca_prior_in$S_n #number of pathogens that tested sensitive - n_obs <- wisca_prior_in$tested_n #number of pathogens tested - # n_tested <- round(n_obs*(wisca_calc_input$prop.sens)) - diff_nr <- n_obs - r #difference between number of pathogens that were tested and those that tested sensitive - - #posterior beta parameters - post_beta_1 <- beta_A + r - post_beta_2 <- beta_B + diff_nr - - priors_df <- bind_cols(gamma_A,B,post_beta_1,post_beta_2) - names(priors_df) <- c("gamma_A","B","post_beta_1","post_beta_2") - return(priors_df) -} - -#simulations per bug -#function to generate simulations -coverage <- function(wisca_sim_in){ - sim_coverage <- wisca_sim_in %>% - mutate(random_incidence = runif(n = nrow(wisca_sim_in), min = 0, max = 1), - random_susceptibility = runif(n = nrow(wisca_sim_in), min = 0, max = 1)) %>% - mutate(simulation_bug = qgamma(random_incidence, shape = gamma_A, scale = B), - simulation_suscep = qbeta(p = random_susceptibility, - shape1 = post_beta_1, shape2 = post_beta_2)) %>% - mutate(sim_inc_total = sum(simulation_bug), - simulated_incidence = simulation_bug/sim_inc_total) %>% - mutate(simulated_bugcover = simulated_incidence *simulation_suscep, - coverage = sum(simulated_bugcover)) - sim_coverage_out <- slice_head(sim_coverage, n = 1)$coverage - return(sim_coverage_out) -} - -##will run this for every regimen -set.seed(1243) -cover <- data.frame() - -for (gg in 1:length(unique(wisca_in$keyantimicrobials))){ - simulation_nr <- 1000 - #generate 1000 coverage estimates and determine mean and CI - wisca <- subset(wisca_in, keyantimicrobials == unique(wisca_in$keyantimicrobials)[gg]) - - if (analysis %in% names(wisca)){ - cover_level <- data.frame() - for (ii in 1:length(levels(wisca[[analysis]]))){ - wisca_prep <- subset(wisca, wisca[[analysis]] == levels(wisca[[analysis]])[ii]) - params_priors <- priors(wisca_prior_in = wisca_prep) - - coverage_simulation_total <- replicate(n = simulation_nr,coverage(wisca_sim_in = params_priors)) - - av_coverage <- mean(coverage_simulation_total) - ci_coverage_lower <- quantile(coverage_simulation_total, probs = 0.025) - ci_coverage_upper <- quantile(coverage_simulation_total, probs = 0.975) - - combine_out <- cbind.data.frame(unique(wisca_in$keyantimicrobials)[gg], - levels(wisca[[analysis]])[ii], - av_coverage, - ci_coverage_lower, - ci_coverage_upper, - row.names = NULL) - names(combine_out) <- c("Regimen",str_to_sentence(analysis),"Coverage","Lower_CI","Upper_CI") - cover_level <- rbind.data.frame(cover_level,combine_out, row.names = NULL) - } - cover <- rbind.data.frame(cover,cover_level, row.names = NULL) - - - } else { - params_priors <- priors(wisca_prior_in = wisca) - - coverage_simulation_total <- replicate(n = simulation_nr,coverage(wisca_sim_in = params_priors)) - - av_coverage <- mean(coverage_simulation_total) - ci_coverage_lower <- quantile(coverage_simulation_total, probs = 0.025) - ci_coverage_upper <- quantile(coverage_simulation_total, probs = 0.975) - - combine_out <- cbind.data.frame(unique(wisca_in$keyantimicrobials)[gg], - av_coverage, - ci_coverage_lower, - ci_coverage_upper, - row.names = NULL) - names(combine_out) <- c("Regimen","Coverage","Lower_CI","Upper_CI") - cover <- rbind.data.frame(cover,combine_out, row.names = NULL) - - } - -} -cover <- cover %>% - mutate(across(.cols = Coverage:Upper_CI, ~ round(.x*100,1))) %>% - rowwise() %>% - mutate(Regimen_full = ifelse(str_detect(Regimen, fixed("+")), - str_c(sapply(unlist(unique(str_split(Regimen, fixed("+")))),ab_name,USE.NAMES = FALSE), collapse = " + "), - ab_name(Regimen))) %>% - relocate(Regimen_full,.before = Regimen) - -return(cover) - -} +# Copyright (c) [2022] [Larisse Bolton (author), Aislinn Cook (contributor)] +# Adapted from: Bielicki JA, Sharland M, Heath PT, et al. Evaluation of the coverage of 3 antibiotic regimens for +# neonatal sepsis in the hospital setting across Asian countries. JAMA Netw Open. +# 2020:3(2):e1921124. doi:10.1001.jamanetworkopen.2019.21124 + + +wisca_funct <- function(wisca_in, analysis) { + # function used to gerenrate coverage estimate parameters from input data + priors <- function(wisca_prior_in) { + # INCIDENCE + # prior Dirichlet (Gamma) parameters + A <- rep(1, times = length(unique(wisca_prior_in$mo))) + Multiplier <- rep(1, times = length(unique(wisca_prior_in$mo))) + multinomial_obs <- round((wisca_prior_in$n.x) * Multiplier, 0) + + # posterior Dirichlet (Gamma) parameters + gamma_A <- A + multinomial_obs + B <- rep(1, times = length(gamma_A)) + + ## SUSCEPTIBILITY + # prior Beta parameters + beta_A <- rep(1, times = length(unique(wisca_prior_in$mo))) + beta_B <- rep(1, times = length(unique(wisca_prior_in$mo))) + + # Binomial distribution parameters + r <- wisca_prior_in$S_n # number of pathogens that tested sensitive + n_obs <- wisca_prior_in$tested_n # number of pathogens tested + # n_tested <- round(n_obs*(wisca_calc_input$prop.sens)) + diff_nr <- n_obs - r # difference between number of pathogens that were tested and those that tested sensitive + + # posterior beta parameters + post_beta_1 <- beta_A + r + post_beta_2 <- beta_B + diff_nr + + priors_df <- bind_cols(gamma_A, B, post_beta_1, post_beta_2) + names(priors_df) <- c("gamma_A", "B", "post_beta_1", "post_beta_2") + return(priors_df) + } + + # simulations per bug + # function to generate simulations + coverage <- function(wisca_sim_in) { + sim_coverage <- wisca_sim_in %>% + mutate( + random_incidence = runif(n = nrow(wisca_sim_in), min = 0, max = 1), + random_susceptibility = runif(n = nrow(wisca_sim_in), min = 0, max = 1) + ) %>% + mutate( + simulation_bug = qgamma(random_incidence, shape = gamma_A, scale = B), + simulation_suscep = qbeta( + p = random_susceptibility, + shape1 = post_beta_1, shape2 = post_beta_2 + ) + ) %>% + mutate( + sim_inc_total = sum(simulation_bug), + simulated_incidence = simulation_bug / sim_inc_total + ) %>% + mutate( + simulated_bugcover = simulated_incidence * simulation_suscep, + coverage = sum(simulated_bugcover) + ) + sim_coverage_out <- slice_head(sim_coverage, n = 1)$coverage + return(sim_coverage_out) + } + + ## will run this for every regimen + set.seed(1243) + cover <- data.frame() + + for (gg in 1:length(unique(wisca_in$keyantimicrobials))) { + simulation_nr <- 1000 + # generate 1000 coverage estimates and determine mean and CI + wisca <- subset(wisca_in, keyantimicrobials == unique(wisca_in$keyantimicrobials)[gg]) + + if (analysis %in% names(wisca)) { + cover_level <- data.frame() + for (ii in 1:length(levels(wisca[[analysis]]))) { + wisca_prep <- subset(wisca, wisca[[analysis]] == levels(wisca[[analysis]])[ii]) + params_priors <- priors(wisca_prior_in = wisca_prep) + + coverage_simulation_total <- replicate(n = simulation_nr, coverage(wisca_sim_in = params_priors)) + + av_coverage <- mean(coverage_simulation_total) + ci_coverage_lower <- quantile(coverage_simulation_total, probs = 0.025) + ci_coverage_upper <- quantile(coverage_simulation_total, probs = 0.975) + + combine_out <- cbind.data.frame(unique(wisca_in$keyantimicrobials)[gg], + levels(wisca[[analysis]])[ii], + av_coverage, + ci_coverage_lower, + ci_coverage_upper, + row.names = NULL + ) + names(combine_out) <- c("Regimen", str_to_sentence(analysis), "Coverage", "Lower_CI", "Upper_CI") + cover_level <- rbind.data.frame(cover_level, combine_out, row.names = NULL) + } + cover <- rbind.data.frame(cover, cover_level, row.names = NULL) + } else { + params_priors <- priors(wisca_prior_in = wisca) + + coverage_simulation_total <- replicate(n = simulation_nr, coverage(wisca_sim_in = params_priors)) + + av_coverage <- mean(coverage_simulation_total) + ci_coverage_lower <- quantile(coverage_simulation_total, probs = 0.025) + ci_coverage_upper <- quantile(coverage_simulation_total, probs = 0.975) + + combine_out <- cbind.data.frame(unique(wisca_in$keyantimicrobials)[gg], + av_coverage, + ci_coverage_lower, + ci_coverage_upper, + row.names = NULL + ) + names(combine_out) <- c("Regimen", "Coverage", "Lower_CI", "Upper_CI") + cover <- rbind.data.frame(cover, combine_out, row.names = NULL) + } + } + cover <- cover %>% + mutate(across(.cols = Coverage:Upper_CI, ~ round(.x * 100, 1))) %>% + rowwise() %>% + mutate(Regimen_full = ifelse(str_detect(Regimen, fixed("+")), + str_c(sapply(unlist(unique(str_split(Regimen, fixed("+")))), ab_name, USE.NAMES = FALSE), collapse = " + "), + ab_name(Regimen) + )) %>% + relocate(Regimen_full, .before = Regimen) + + return(cover) +} diff --git a/index.md b/index.md index 4f4ab3797..401da2dd4 100644 --- a/index.md +++ b/index.md @@ -170,13 +170,9 @@ example_isolates %>% #> ℹ Using column mo as input for `mo_fullname()` #> ℹ Using column mo as input for `mo_is_gram_negative()` #> ℹ Using column mo as input for `mo_is_intrinsic_resistant()` -#> ℹ Determining intrinsic resistance based on 'EUCAST Expected -#> Resistant Phenotypes' v1.2 (2023). This note will be shown -#> once per session. -#> ℹ For `aminoglycosides()` using columns GEN (gentamicin), TOB -#> (tobramycin), AMK (amikacin), and KAN (kanamycin) -#> ℹ For `carbapenems()` using columns IPM (imipenem) and MEM -#> (meropenem) +#> ℹ Determining intrinsic resistance based on 'EUCAST Expected Resistant Phenotypes' v1.2 (2023). This note will be shown once per session. +#> ℹ For `aminoglycosides()` using columns GEN (gentamicin), TOB (tobramycin), AMK (amikacin), and KAN (kanamycin) +#> ℹ For `carbapenems()` using columns IPM (imipenem) and MEM (meropenem) #> # A tibble: 35 × 7 #> bacteria GEN TOB AMK KAN IPM MEM #> @@ -253,8 +249,7 @@ for tracking resistance per species over time: antibiogram(example_isolates, mo_transform = "gramstain", antimicrobials = c("AMC", carbapenems(), "TZP")) -#> ℹ For `carbapenems()` using columns IPM (imipenem) and MEM -#> (meropenem) +#> ℹ For `carbapenems()` using columns IPM (imipenem) and MEM (meropenem) ``` | Pathogen | Amoxicillin/clavulanic acid | Imipenem | Meropenem | Piperacillin/tazobactam | @@ -369,16 +364,13 @@ out <- example_isolates %>% # calculate AMR using resistance(), over all aminoglycosides and polymyxins: summarise(across(c(aminoglycosides(), polymyxins()), resistance)) -#> ℹ For `aminoglycosides()` using columns GEN (gentamicin), TOB -#> (tobramycin), AMK (amikacin), and KAN (kanamycin) +#> ℹ For `aminoglycosides()` using columns GEN (gentamicin), TOB (tobramycin), AMK (amikacin), and KAN (kanamycin) #> ℹ For `polymyxins()` using column COL (colistin) #> Warning: There was 1 warning in `summarise()`. -#> ℹ In argument: `across(c(aminoglycosides(), polymyxins()), -#> resistance)`. +#> ℹ In argument: `across(c(aminoglycosides(), polymyxins()), resistance)`. #> ℹ In group 3: `ward = "Outpatient"`. #> Caused by warning: -#> ! Introducing NA: only 23 results available for KAN in group: -#> ward = "Outpatient" (whilst `minimum = 30`). +#> ! Introducing NA: only 23 results available for KAN in group: ward = "Outpatient" (whilst `minimum = 30`). out #> # A tibble: 3 × 6 #> ward GEN TOB AMK KAN COL diff --git a/man/antibiogram.Rd b/man/antibiogram.Rd index 1124b92e7..8e563d409 100644 --- a/man/antibiogram.Rd +++ b/man/antibiogram.Rd @@ -334,7 +334,8 @@ if (requireNamespace("dplyr")) { top_n_microorganisms(n = 10) \%>\% group_by( age_group = age_groups(age, c(25, 50, 75)), - gender) \%>\% + gender + ) \%>\% wisca(antimicrobials = c("TZP", "TZP+TOB", "TZP+GEN")) } diff --git a/man/antimicrobial_selectors.Rd b/man/antimicrobial_selectors.Rd index 5b9005a67..359c9ba47 100644 --- a/man/antimicrobial_selectors.Rd +++ b/man/antimicrobial_selectors.Rd @@ -308,7 +308,6 @@ carbapenems() # Though they are primarily intended to use for selections and filters. # Examples sections below are split into 'dplyr', 'base R', and 'data.table': -\donttest{ \dontrun{ # dplyr ------------------------------------------------------------------- @@ -329,8 +328,10 @@ example_isolates \%>\% select(administrable_per_os()) # get AMR for all aminoglycosides e.g., per ward: example_isolates \%>\% group_by(ward) \%>\% - summarise(across(aminoglycosides(), - resistance)) + summarise(across( + aminoglycosides(), + resistance + )) # You can combine selectors with '&' to be more specific: example_isolates \%>\% @@ -340,8 +341,10 @@ example_isolates \%>\% example_isolates \%>\% filter(mo_genus() \%in\% c("Escherichia", "Klebsiella")) \%>\% group_by(ward) \%>\% - summarise_at(not_intrinsic_resistant(), - resistance) + summarise_at( + not_intrinsic_resistant(), + resistance + ) # get susceptibility for antimicrobials whose name contains "trim": example_isolates \%>\% @@ -394,8 +397,9 @@ x <- example_isolates[carbapenems() == "R", ] y <- example_isolates \%>\% filter(carbapenems() == "R") z <- example_isolates \%>\% filter(if_all(carbapenems(), ~ .x == "R")) identical(x, y) && identical(y, z) - } + +\donttest{ # base R ------------------------------------------------------------------ # select columns 'IPM' (imipenem) and 'MEM' (meropenem)