From 046d19506469559f4b3d658387679496e38cede6 Mon Sep 17 00:00:00 2001 From: "Matthijs S. Berends" Date: Tue, 15 Jan 2019 12:45:24 +0100 Subject: [PATCH] resistance predict --- DESCRIPTION | 2 +- NAMESPACE | 4 + NEWS.md | 6 + R/eucast_rules.R | 5 +- R/first_isolate.R | 27 +-- R/key_antibiotics.R | 6 +- R/mdro.R | 5 +- R/misc.R | 45 ++++ R/resistance_predict.R | 261 ++++++++++++++--------- docs/news/index.html | 2 + docs/reference/first_isolate.html | 2 +- docs/reference/index.html | 2 +- docs/reference/resistance_predict.html | 71 +++--- man/first_isolate.Rd | 2 +- man/resistance_predict.Rd | 70 +++--- tests/testthat/test-resistance_predict.R | 5 + 16 files changed, 306 insertions(+), 209 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 9a86ca29..aaea9a0d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: AMR Version: 0.5.0.9009 -Date: 2019-01-12 +Date: 2019-01-15 Title: Antimicrobial Resistance Analysis Authors@R: c( person( diff --git a/NAMESPACE b/NAMESPACE index 3ce45c98..e5ea1a50 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -20,6 +20,7 @@ S3method(kurtosis,default) S3method(kurtosis,matrix) S3method(plot,frequency_tbl) S3method(plot,mic) +S3method(plot,resistance_predict) S3method(plot,rsi) S3method(print,atc) S3method(print,frequency_tbl) @@ -76,6 +77,7 @@ export(g.test) export(geom_rsi) export(get_locale) export(ggplot_rsi) +export(ggplot_rsi_predict) export(guess_ab_col) export(guess_atc) export(guess_mo) @@ -215,6 +217,7 @@ importFrom(dplyr,mutate_all) importFrom(dplyr,mutate_at) importFrom(dplyr,n) importFrom(dplyr,n_distinct) +importFrom(dplyr,n_groups) importFrom(dplyr,progress_estimated) importFrom(dplyr,pull) importFrom(dplyr,row_number) @@ -228,6 +231,7 @@ importFrom(dplyr,top_n) importFrom(dplyr,ungroup) importFrom(dplyr,vars) importFrom(grDevices,boxplot.stats) +importFrom(graphics,arrows) importFrom(graphics,axis) importFrom(graphics,barplot) importFrom(graphics,hist) diff --git a/NEWS.md b/NEWS.md index d40f259b..fa50d777 100755 --- a/NEWS.md +++ b/NEWS.md @@ -12,6 +12,12 @@ * Function `mo_renamed()` to get a list of all returned values from `as.mo()` that have had taxonomic renaming * Function `age()` to calculate the (patients) age in years * Function `age_groups()` to split ages into custom or predefined groups (like children or elderly). This allows for easier demographic antimicrobial resistance analysis per age group. +* Function `ggplot_rsi_predict()` as well as the base R `plot()` function can now be used for resistance prediction calculated with `resistance_predict()`: + ```r + x <- resistance_predict(septic_patients, col_ab = "amox") + plot(x) + ggplot_rsi_predict(x) + ``` * Functions `filter_first_isolate()` and `filter_first_weighted_isolate()` to shorten and fasten filtering on data sets with antimicrobial results, e.g.: ```r septic_patients %>% filter_first_isolate(...) diff --git a/R/eucast_rules.R b/R/eucast_rules.R index 6ae03d7d..459db3cf 100755 --- a/R/eucast_rules.R +++ b/R/eucast_rules.R @@ -231,9 +231,8 @@ eucast_rules <- function(tbl, # try to find columns based on type # -- mo - if (is.null(col_mo) & "mo" %in% lapply(tbl, class)) { - col_mo <- colnames(tbl)[lapply(tbl, class) == "mo"][1] - message(blue(paste0("NOTE: Using column `", bold(col_mo), "` as input for `col_mo`."))) + if (is.null(col_mo)) { + col_mo <- search_type_in_df(tbl = tbl, type = "mo") } if (is.null(col_mo)) { stop("`col_mo` must be set.", call. = FALSE) diff --git a/R/first_isolate.R b/R/first_isolate.R index 2694caf7..20335186 100755 --- a/R/first_isolate.R +++ b/R/first_isolate.R @@ -23,7 +23,7 @@ #' #' Determine first (weighted) isolates of all microorganisms of every patient per episode and (if needed) per specimen type. #' @param tbl a \code{data.frame} containing isolates. -#' @param col_date column name of the result date (or date that is was received on the lab), defaults to the first column of class \code{Date} +#' @param col_date column name of the result date (or date that is was received on the lab), defaults to the first column of with a date class #' @param col_patient_id column name of the unique IDs of the patients, defaults to the first column that starts with 'patient' or 'patid' (case insensitive) #' @param col_mo column name of the unique IDs of the microorganisms (see \code{\link{mo}}), defaults to the first column of class \code{mo}. Values will be coerced using \code{\link{as.mo}}. #' @param col_testcode column name of the test codes. Use \code{col_testcode = NULL} to \strong{not} exclude certain test codes (like test codes for screening). In that case \code{testcodes_exclude} will be ignored. @@ -187,9 +187,8 @@ first_isolate <- function(tbl, # try to find columns based on type # -- mo - if (is.null(col_mo) & "mo" %in% lapply(tbl, class)) { - col_mo <- colnames(tbl)[lapply(tbl, class) == "mo"][1] - message(blue(paste0("NOTE: Using column `", bold(col_mo), "` as input for `col_mo`."))) + if (is.null(col_mo)) { + col_mo <- search_type_in_df(tbl = tbl, type = "mo") } if (is.null(col_mo)) { stop("`col_mo` must be set.", call. = FALSE) @@ -197,33 +196,25 @@ first_isolate <- function(tbl, # -- date if (is.null(col_date)) { - for (i in 1:ncol(tbl)) { - if ("Date" %in% class(tbl %>% pull(i)) | "POSIXct" %in% class(tbl %>% pull(i))) { - col_date <- colnames(tbl)[i] - message(blue(paste0("NOTE: Using column `", bold(col_date), "` as input for `col_date`."))) - break - } - } + col_date <- search_type_in_df(tbl = tbl, type = "date") } if (is.null(col_date)) { stop("`col_date` must be set.", call. = FALSE) } - # convert to Date (pipes for supporting tibbles too) + # convert to Date (pipes/pull for supporting tibbles too) tbl[, col_date] <- tbl %>% pull(col_date) %>% as.Date() # -- patient id - if (is.null(col_patient_id) & any(colnames(tbl) %like% "^(patient|patid)")) { - col_patient_id <- colnames(tbl)[colnames(tbl) %like% "^(patient|patid)"][1] - message(blue(paste0("NOTE: Using column `", bold(col_patient_id), "` as input for `col_patient_id`."))) + if (is.null(col_patient_id)) { + col_patient_id <- search_type_in_df(tbl = tbl, type = "patient_id") } if (is.null(col_patient_id)) { stop("`col_patient_id` must be set.", call. = FALSE) } # -- key antibiotics - if (is.null(col_keyantibiotics) & any(colnames(tbl) %like% "^key.*(ab|antibiotics)")) { - col_keyantibiotics <- colnames(tbl)[colnames(tbl) %like% "^key.*(ab|antibiotics)"][1] - message(blue(paste0("NOTE: Using column `", bold(col_keyantibiotics), "` as input for `col_keyantibiotics`. Use ", bold("col_keyantibiotics = FALSE"), " to prevent this."))) + if (is.null(col_keyantibiotics)) { + col_keyantibiotics <- search_type_in_df(tbl = tbl, type = "keyantibiotics") } if (isFALSE(col_keyantibiotics)) { col_keyantibiotics <- NULL diff --git a/R/key_antibiotics.R b/R/key_antibiotics.R index f654c01f..b5804a11 100755 --- a/R/key_antibiotics.R +++ b/R/key_antibiotics.R @@ -101,9 +101,8 @@ key_antibiotics <- function(tbl, # try to find columns based on type # -- mo - if (is.null(col_mo) & "mo" %in% lapply(tbl, class)) { - col_mo <- colnames(tbl)[lapply(tbl, class) == "mo"][1] - message(blue(paste0("NOTE: Using column `", bold(col_mo), "` as input for `col_mo`."))) + if (is.null(col_mo)) { + col_mo <- search_type_in_df(tbl = tbl, type = "mo") } if (is.null(col_mo)) { stop("`col_mo` must be set.", call. = FALSE) @@ -114,7 +113,6 @@ key_antibiotics <- function(tbl, GramPos_1, GramPos_2, GramPos_3, GramPos_4, GramPos_5, GramPos_6, GramNeg_1, GramNeg_2, GramNeg_3, GramNeg_4, GramNeg_5, GramNeg_6) col.list <- check_available_columns(tbl = tbl, col.list = col.list, info = warnings) - print(col.list) universal_1 <- col.list[universal_1] universal_2 <- col.list[universal_2] universal_3 <- col.list[universal_3] diff --git a/R/mdro.R b/R/mdro.R index ae139bc8..51a9c0f8 100755 --- a/R/mdro.R +++ b/R/mdro.R @@ -113,9 +113,8 @@ mdro <- function(tbl, # try to find columns based on type # -- mo - if (is.null(col_mo) & "mo" %in% lapply(tbl, class)) { - col_mo <- colnames(tbl)[lapply(tbl, class) == "mo"][1] - message(blue(paste0("NOTE: Using column `", bold(col_mo), "` as input for `col_mo`."))) + if (is.null(col_mo)) { + col_mo <- search_type_in_df(tbl = tbl, type = "mo") } if (is.null(col_mo)) { stop("`col_mo` must be set.", call. = FALSE) diff --git a/R/misc.R b/R/misc.R index c2b17c6a..846ebefe 100755 --- a/R/misc.R +++ b/R/misc.R @@ -123,3 +123,48 @@ size_humanreadable <- function(bytes, decimals = 1) { out <- paste(sprintf(paste0("%.", decimals, "f"), bytes / (1024 ^ factor)), size[factor + 1]) out } + +#' @importFrom crayon blue bold +#' @importFrom dplyr %>% pull +search_type_in_df <- function(tbl, type) { + # try to find columns based on type + found <- NULL + + # -- mo + if (type == "mo") { + if ("mo" %in% lapply(tbl, class)) { + found <- colnames(tbl)[lapply(tbl, class) == "mo"][1] + } + } + # -- key antibiotics + if (type == "keyantibiotics") { + if (any(colnames(tbl) %like% "^key.*(ab|antibiotics)")) { + found <- colnames(tbl)[colnames(tbl) %like% "^key.*(ab|antibiotics)"][1] + } + } + # -- date + if (type == "date") { + for (i in 1:ncol(tbl)) { + if ("Date" %in% class(tbl %>% pull(i)) | "POSIXct" %in% class(tbl %>% pull(i))) { + found <- colnames(tbl)[i] + break + } + } + + } + # -- patient id + if (type == "patient_id") { + if (any(colnames(tbl) %like% "^(patient|patid)")) { + found <- colnames(tbl)[colnames(tbl) %like% "^(patient|patid)"][1] + } + } + + if (!is.null(found)) { + msg <- paste0("NOTE: Using column `", bold(found), "` as input for `col_", type, "`.") + if (type == "keyantibiotics") { + msg <- paste(msg, "Use", bold("col_keyantibiotics = FALSE"), "to prevent this.") + } + message(blue(msg)) + } + found +} diff --git a/R/resistance_predict.R b/R/resistance_predict.R index 6caab1d0..bf2ea9e7 100755 --- a/R/resistance_predict.R +++ b/R/resistance_predict.R @@ -23,17 +23,18 @@ #' #' Create a prediction model to predict antimicrobial resistance for the next years on statistical solid ground. Standard errors (SE) will be returned as columns \code{se_min} and \code{se_max}. See Examples for a real live example. #' @inheritParams first_isolate +#' @inheritParams graphics::plot #' @param col_ab column name of \code{tbl} with antimicrobial interpretations (\code{R}, \code{I} and \code{S}) -#' @param col_date column name of the date, will be used to calculate years if this column doesn't consist of years already +#' @param col_date column name of the date, will be used to calculate years if this column doesn't consist of years already, defaults to the first column of with a date class #' @param year_min lowest year to use in the prediction model, dafaults to the lowest year in \code{col_date} #' @param year_max highest year to use in the prediction model, defaults to 10 years after today #' @param year_every unit of sequence between lowest year found in the data and \code{year_max} #' @param minimum minimal amount of available isolates per year to include. Years containing less observations will be estimated by the model. -#' @param model the statistical model of choice. Valid values are \code{"binomial"} (or \code{"binom"} or \code{"logit"}) or \code{"loglin"} or \code{"linear"} (or \code{"lin"}). +#' @param model the statistical model of choice. Valid values are \code{"binomial"} (or \code{"binom"} or \code{"logit"}) or \code{"loglin"} (or \code{"poisson"}) or \code{"linear"} (or \code{"lin"}). #' @param I_as_R a logical to indicate whether values \code{I} should be treated as \code{R} #' @param preserve_measurements a logical to indicate whether predictions of years that are actually available in the data should be overwritten by the original data. The standard errors of those years will be \code{NA}. #' @param info a logical to indicate whether textual analysis should be printed with the name and \code{\link{summary}} of the statistical model. -#' @return \code{data.frame} with columns: +#' @return \code{data.frame} with extra class \code{"resistance_predict"} with columns: #' \itemize{ #' \item{\code{year}} #' \item{\code{value}, the same as \code{estimated} when \code{preserve_measurements = FALSE}, and a combination of \code{observed} and \code{estimated} otherwise} @@ -47,42 +48,23 @@ #' @rdname resistance_predict #' @export #' @importFrom stats predict glm lm -#' @importFrom dplyr %>% pull mutate mutate_at n group_by_at summarise filter filter_at all_vars n_distinct arrange case_when +#' @importFrom dplyr %>% pull mutate mutate_at n group_by_at summarise filter filter_at all_vars n_distinct arrange case_when n_groups #' @inheritSection AMR Read more on our website! #' @examples -#' \dontrun{ -#' # use it with base R: -#' resistance_predict(tbl = tbl[which(first_isolate == TRUE & genus == "Haemophilus"),], -#' col_ab = "amcl", col_date = "date") +#' x <- resistance_predict(septic_patients, col_ab = "amox", year_min = 2010) +#' plot(x) +#' ggplot_rsi_predict(x) #' -#' # or use dplyr so you can actually read it: +#' # use dplyr so you can actually read it: #' library(dplyr) -#' tbl %>% -#' filter(first_isolate == TRUE, -#' genus == "Haemophilus") %>% -#' resistance_predict(amcl, date) -#' } +#' x <- septic_patients %>% +#' filter_first_isolate() %>% +#' filter(mo_genus(mo) == "Staphylococcus") %>% +#' resistance_predict("peni") +#' plot(x) #' #' -#' # real live example: -#' library(dplyr) -#' septic_patients %>% -#' # get bacteria properties like genus and species -#' left_join_microorganisms("mo") %>% -#' # calculate first isolates -#' mutate(first_isolate = first_isolate(.)) %>% -#' # filter on first E. coli isolates -#' filter(genus == "Escherichia", -#' species == "coli", -#' first_isolate == TRUE) %>% -#' # predict resistance of cefotaxime for next years -#' resistance_predict(col_ab = "cfot", -#' col_date = "date", -#' year_max = 2025, -#' preserve_measurements = TRUE, -#' minimum = 0) -#' -#' # create nice plots with ggplot +#' # create nice plots with ggplot yourself #' if (!require(ggplot2)) { #' #' data <- septic_patients %>% @@ -90,7 +72,7 @@ #' resistance_predict(col_ab = "amox", #' col_date = "date", #' info = FALSE, -#' minimum = 15) +#' minimum = 15) #' #' ggplot(data, #' aes(x = year)) + @@ -110,7 +92,7 @@ #' } resistance_predict <- function(tbl, col_ab, - col_date, + col_date = NULL, year_min = NULL, year_max = NULL, year_every = 1, @@ -128,23 +110,23 @@ resistance_predict <- function(tbl, stop('Column ', col_ab, ' not found.') } + # -- date + if (is.null(col_date)) { + col_date <- search_type_in_df(tbl = tbl, type = "date") + } + if (is.null(col_date)) { + stop("`col_date` must be set.", call. = FALSE) + } + if (!col_date %in% colnames(tbl)) { stop('Column ', col_date, ' not found.') } - if ('grouped_df' %in% class(tbl)) { + + if (n_groups(tbl) > 1) { # no grouped tibbles please, mutate will throw errors tbl <- base::as.data.frame(tbl, stringsAsFactors = FALSE) } - if (I_as_R == TRUE) { - tbl[, col_ab] <- gsub('I', 'R', tbl %>% pull(col_ab)) - } - - tbl <- tbl %>% - mutate_at(col_ab, as.rsi) %>% - filter_at(col_ab, all_vars(!is.na(.))) - tbl[, col_ab] <- droplevels(tbl[, col_ab]) - year <- function(x) { if (all(grepl('^[0-9]{4}$', x))) { x @@ -154,13 +136,23 @@ resistance_predict <- function(tbl, } df <- tbl %>% - mutate(year = tbl %>% pull(col_date) %>% year()) %>% + mutate_at(col_ab, as.rsi) %>% + mutate_at(col_ab, droplevels) %>% + mutate_at(col_ab, funs( + if (I_as_R == TRUE) { + gsub("I", "R", .) + } else { + gsub("I", "S", .) + } + )) %>% + filter_at(col_ab, all_vars(!is.na(.))) %>% + mutate(year = pull(., col_date) %>% year()) %>% group_by_at(c('year', col_ab)) %>% summarise(n()) if (df %>% pull(col_ab) %>% n_distinct(na.rm = TRUE) < 2) { stop("No variety in antimicrobial interpretations - all isolates are '", - df %>% pull(col_ab) %>% unique() %>% .[!is.na(.)], "'.", + df %>% pull(col_ab) %>% unique(), "'.", call. = FALSE) } @@ -168,8 +160,11 @@ resistance_predict <- function(tbl, df <- df %>% filter(!is.na(antibiotic)) %>% tidyr::spread(antibiotic, observations, fill = 0) %>% - mutate(total = R + S) %>% - filter(total >= minimum) + filter((R + S) >= minimum) + df_matrix <- df %>% + ungroup() %>% + select(R, S) %>% + as.matrix() if (NROW(df) == 0) { stop('There are no observations.') @@ -185,41 +180,44 @@ resistance_predict <- function(tbl, year_max <- year(Sys.Date()) + 10 } - years_predict <- seq(from = year_min, to = year_max, by = year_every) + years <- list(year = seq(from = year_min, to = year_max, by = year_every)) if (model %in% c('binomial', 'binom', 'logit')) { - logitmodel <- with(df, glm(cbind(R, S) ~ year, family = binomial)) + model <- "binomial" + model_lm <- with(df, glm(df_matrix ~ year, family = binomial)) if (info == TRUE) { cat('\nLogistic regression model (logit) with binomial distribution') cat('\n------------------------------------------------------------\n') - print(summary(logitmodel)) + print(summary(model_lm)) } - predictmodel <- predict(logitmodel, newdata = with(df, list(year = years_predict)), type = "response", se.fit = TRUE) + predictmodel <- predict(model_lm, newdata = years, type = "response", se.fit = TRUE) prediction <- predictmodel$fit se <- predictmodel$se.fit - } else if (model == 'loglin') { - loglinmodel <- with(df, glm(R ~ year, family = poisson)) + } else if (model %in% c('loglin', 'poisson')) { + model <- "poisson" + model_lm <- with(df, glm(R ~ year, family = poisson)) if (info == TRUE) { cat('\nLog-linear regression model (loglin) with poisson distribution') cat('\n--------------------------------------------------------------\n') - print(summary(loglinmodel)) + print(summary(model_lm)) } - predictmodel <- predict(loglinmodel, newdata = with(df, list(year = years_predict)), type = "response", se.fit = TRUE) + predictmodel <- predict(model_lm, newdata = years, type = "response", se.fit = TRUE) prediction <- predictmodel$fit se <- predictmodel$se.fit } else if (model %in% c('lin', 'linear')) { - linmodel <- with(df, lm((R / (R + S)) ~ year)) + model <- "linear" + model_lm <- with(df, lm((R / (R + S)) ~ year)) if (info == TRUE) { cat('\nLinear regression model') cat('\n-----------------------\n') - print(summary(linmodel)) + print(summary(model_lm)) } - predictmodel <- predict(linmodel, newdata = with(df, list(year = years_predict)), se.fit = TRUE) + predictmodel <- predict(model_lm, newdata = years, se.fit = TRUE) prediction <- predictmodel$fit se <- predictmodel$se.fit @@ -228,64 +226,117 @@ resistance_predict <- function(tbl, } # prepare the output dataframe - prediction <- data.frame(year = years_predict, value = prediction, stringsAsFactors = FALSE) + df_prediction <- data.frame(year = unlist(years), + value = prediction, + stringsAsFactors = FALSE) %>% - prediction$se_min <- prediction$value - se - prediction$se_max <- prediction$value + se + mutate(se_min = value - se, + se_max = value + se) - if (model == 'loglin') { - prediction$value <- prediction$value %>% - format(scientific = FALSE) %>% - as.integer() - prediction$se_min <- prediction$se_min %>% as.integer() - prediction$se_max <- prediction$se_max %>% as.integer() - - colnames(prediction) <- c('year', 'amountR', 'se_max', 'se_min') + if (model == 'poisson') { + df_prediction <- df_prediction %>% + mutate(value = value %>% + format(scientific = FALSE) %>% + as.integer(), + se_min = as.integer(se_min), + se_max = as.integer(se_max)) } else { - prediction$se_max[which(prediction$se_max > 1)] <- 1 + df_prediction <- df_prediction %>% + # se_max not above 1 + mutate(se_max = ifelse(se_max > 1, 1, se_max)) } - prediction$se_min[which(prediction$se_min < 0)] <- 0 - prediction$observations = NA + df_prediction <- df_prediction %>% + # se_min not below 0 + mutate(se_min = ifelse(se_min < 0, 0, se_min)) - total <- prediction + df_observations <- df %>% + ungroup() %>% + transmute(year, + observations = R + S, + observed = R / (R + S)) + df_prediction <- df_prediction %>% + left_join(df_observations, by = "year") %>% + mutate(estimated = value) if (preserve_measurements == TRUE) { # replace estimated data by observed data - if (I_as_R == TRUE) { - if (!'I' %in% colnames(df)) { - df$I <- 0 - } - df$value <- df$R / rowSums(df[, c('R', 'S', 'I')]) - } else { - df$value <- df$R / rowSums(df[, c('R', 'S')]) - } - measurements <- data.frame(year = df$year, - value = df$value, - se_min = NA, - se_max = NA, - observations = df$total, - stringsAsFactors = FALSE) - colnames(measurements) <- colnames(prediction) - - total <- rbind(measurements, - prediction %>% filter(!year %in% df$year)) - if (model %in% c('binomial', 'binom', 'logit')) { - total <- total %>% mutate(observed = ifelse(is.na(observations), NA, value), - estimated = prediction$value) - } + df_prediction <- df_prediction %>% + mutate(value = ifelse(!is.na(observed), observed, value), + se_min = ifelse(!is.na(observed), NA, se_min), + se_max = ifelse(!is.na(observed), NA, se_max)) } - if ("value" %in% colnames(total)) { - total <- total %>% - mutate(value = case_when(value > 1 ~ 1, - value < 0 ~ 0, - TRUE ~ value)) - } - - total %>% arrange(year) + df_prediction <- df_prediction %>% + mutate(value = case_when(value > 1 ~ 1, + value < 0 ~ 0, + TRUE ~ value)) %>% + arrange(year) + structure( + .Data = df_prediction, + class = c("resistance_predict", "data.frame"), + I_as_R = I_as_R, + model_title = model, + model = model_lm, + ab = col_ab + ) } #' @rdname resistance_predict #' @export rsi_predict <- resistance_predict + +#' @exportMethod plot.mic +#' @export +#' @importFrom dplyr %>% group_by summarise +#' @importFrom graphics plot axis arrows +#' @rdname resistance_predict +plot.resistance_predict <- function(x, main = paste("Resistance prediction of", attributes(x)$ab), ...) { + if (attributes(x)$I_as_R == TRUE) { + ylab <- "%IR" + } else { + ylab <- "%R" + } + plot(x = x$year, + y = x$value, + ylim = c(0, 1), + yaxt = "n", # no y labels + pch = 19, # closed dots + ylab = paste0("Percentage (", ylab, ")"), + xlab = "Year", + main = main, + sub = paste0("(model: ", attributes(x)$model_title, ")")) + + axis(side = 2, at = seq(0, 1, 0.1), labels = paste0(0:10 * 10, "%")) + + # arrows hack: https://stackoverflow.com/a/22037078/4575331 + arrows(x0 = x$year, + y0 = x$se_min, + x1 = x$year, + y1 = x$se_max, length = 0.05, angle = 90, code = 3) +} + +#' @rdname resistance_predict +#' @export +ggplot_rsi_predict <- function(x, main = paste("Resistance prediction of", attributes(x)$ab), ...) { + + if (!"resistance_predict" %in% class(x)) { + stop("`x` must be a resistance prediction model created with resistance_predict().") + } + + if (attributes(x)$I_as_R == TRUE) { + ylab <- "%IR" + } else { + ylab <- "%R" + } + suppressWarnings( + ggplot2::ggplot(x, ggplot2::aes(x = year, y = value)) + + ggplot2::geom_col() + + ggplot2::geom_errorbar(ggplot2::aes(ymin = se_min, ymax = se_max)) + + scale_y_percent() + + labs(title = main, + y = paste0("Percentage (", ylab, ")"), + x = "Year", + caption = paste0("(model: ", attributes(x)$model_title, ")")) + ) +} diff --git a/docs/news/index.html b/docs/news/index.html index 25866add..e1d1981c 100644 --- a/docs/news/index.html +++ b/docs/news/index.html @@ -238,6 +238,8 @@
  • Function mo_renamed() to get a list of all returned values from as.mo() that have had taxonomic renaming
  • Function age() to calculate the (patients) age in years
  • Function age_groups() to split ages into custom or predefined groups (like children or elderly). This allows for easier demographic antimicrobial resistance analysis per age group.
    • +
  • Function ggplot_rsi_predict() as well as the base R plot() function can now be used for resistance prediction calculated with resistance_predict(): r x <- resistance_predict(septic_patients, col_ab = "amox") plot(x) ggplot_rsi_predict(x) +
  • Functions filter_first_isolate() and filter_first_weighted_isolate() to shorten and fasten filtering on data sets with antimicrobial results, e.g.: r septic_patients %>% filter_first_isolate(...) # or filter_first_isolate(septic_patients, ...) is equal to: r septic_patients %>% mutate(only_firsts = first_isolate(septic_patients, ...)) %>% filter(only_firsts == TRUE) %>% select(-only_firsts)
  • New vignettes about how to conduct AMR analysis, predict antimicrobial resistance, use the G-test and more. These are also available (and even easier readable) on our website: https://msberends.gitlab.io/AMR.
    • diff --git a/docs/reference/first_isolate.html b/docs/reference/first_isolate.html index a45519bf..e100c72a 100644 --- a/docs/reference/first_isolate.html +++ b/docs/reference/first_isolate.html @@ -250,7 +250,7 @@ col_date -

    column name of the result date (or date that is was received on the lab), defaults to the first column of class Date

    +

    column name of the result date (or date that is was received on the lab), defaults to the first column of with a date class

    col_patient_id diff --git a/docs/reference/index.html b/docs/reference/index.html index a85507b4..3b8334ee 100644 --- a/docs/reference/index.html +++ b/docs/reference/index.html @@ -415,7 +415,7 @@ -

    resistance_predict() rsi_predict()

    +

    resistance_predict() rsi_predict() plot(<resistance_predict>) ggplot_rsi_predict()

    Predict antimicrobial resistance

    diff --git a/docs/reference/resistance_predict.html b/docs/reference/resistance_predict.html index e727669f..4a9d5ab7 100644 --- a/docs/reference/resistance_predict.html +++ b/docs/reference/resistance_predict.html @@ -227,14 +227,22 @@ - 
    resistance_predict(tbl, col_ab, col_date, year_min = NULL,
+    
    resistance_predict(tbl, col_ab, col_date = NULL, year_min = NULL,
   year_max = NULL, year_every = 1, minimum = 30,
   model = "binomial", I_as_R = TRUE, preserve_measurements = TRUE,
   info = TRUE)
 
-rsi_predict(tbl, col_ab, col_date, year_min = NULL, year_max = NULL,
-  year_every = 1, minimum = 30, model = "binomial", I_as_R = TRUE,
-  preserve_measurements = TRUE, info = TRUE)
    
+rsi_predict(tbl, col_ab, col_date = NULL, year_min = NULL,
+  year_max = NULL, year_every = 1, minimum = 30,
+  model = "binomial", I_as_R = TRUE, preserve_measurements = TRUE,
+  info = TRUE)
+
+# S3 method for resistance_predict
+plot(x,
+  main = paste("Resistance prediction of", attributes(x)$ab), ...)
+
+ggplot_rsi_predict(x, main = paste("Resistance prediction of",
+  attributes(x)$ab), ...)

    Arguments

    @@ -249,7 +257,7 @@ - + @@ -269,7 +277,7 @@ - + @@ -283,11 +291,21 @@ + + + + + + + +
    col_date

    column name of the date, will be used to calculate years if this column doesn't consist of years already

    column name of the date, will be used to calculate years if this column doesn't consist of years already, defaults to the first column of with a date class
    year_min
    model

    the statistical model of choice. Valid values are "binomial" (or "binom" or "logit") or "loglin" or "linear" (or "lin").

    the statistical model of choice. Valid values are "binomial" (or "binom" or "logit") or "loglin" (or "poisson") or "linear" (or "lin").
    I_as_Rinfo

    a logical to indicate whether textual analysis should be printed with the name and summary of the statistical model.
    x

    the coordinates of points in the plot. Alternatively, a + single plotting structure, function or any R object with a + plot method can be provided.

    ...

    parameters passed on to the first_isolate function

    Value

    -

    data.frame with columns: