R/plot_readmissions_timeseries.R

In HorridTom/hospitalflow: Analysis of Patient Flow in Hospitals

#' Readmission frequency over period of time
#'
#' @description
#' \code{plot_readmissions_timeseries} shows readmission percentage over period
#' of time. Indicates violation of Schewart rules as well.
#'
#' @param data Hospital episode data.
#' @param startDate The first date of the period for which the analysis will run.
#' @param endDate The last date of the period for which the analysis will run.
#' @param readmissionBy Integer indicating readmission timeframe.
#' @param returnPlot Plots the graph if set to TRUE, returns a dataframe otherwise.
#' @param hospitalName The name of the hospital for which the analysis is being
#' done.
#'
#' @return A lineplot (default) or a dataframe showing the proportion of readmissions
#' over provided period of time.
#'
#' @examples
#' \dontrun{
#' TBD
#' }
#' @export
plot_readmissions_timeseries <- function(data,
                                         startDate,
                                         endDate,
                                         readmissionBy,
                                         returnPlot = TRUE,
                                         hospitalName = "{hospital_name}") {

  # get time zone of data
  time_zone <- attr(data$spell_start, "tzone")

  # set input dates to have the same time zone as the data
  start_date <- as.POSIXct(startDate, tz = time_zone)
  end_date <- as.POSIXct(endDate, tz = time_zone)

  dt_select <- data %>%
    dplyr::select(pseudo_id, spell_number, spell_start, spell_end, ed_admission, admission_method_type) %>%
    dplyr::filter(start_date <= spell_end & end_date >= spell_start) %>%
    dplyr::filter(admission_method_type == "Emergency Admissions" & ed_admission == TRUE)

  #################################################################################################################
  # first, we create a table that contains the csn of the relevant cases
  # and the days since last admission

  dt_calc <- dt_select %>% # we take our data frame
    dplyr::group_by(pseudo_id) %>%
    dplyr::arrange(spell_start) %>%
    dplyr::mutate(
      readm_date = dplyr::lead(spell_start),
      time_to_readm = difftime(readm_date, spell_end, units = c("days")),
      readmission = time_to_readm <= readmissionBy
    ) %>%
    dplyr::ungroup() %>%
    dplyr::filter(spell_end < lubridate::floor_date(
      end_date - lubridate::duration(readmissionBy, units = "days"),
      unit = "month"
    )) %>%
    dplyr::mutate(one_month = lubridate::floor_date(spell_end, "1 month", "month"))

  if (nrow(dt_calc) == 0) {
    warning("Insufficient data to calculate monthly readmission rate.
            Provide more than one month's worth of spell data.")
    return(NULL)
  }

  dt_calc_disch <- dt_calc %>%
    dplyr::group_by(one_month) %>%
    dplyr::summarise(N = dplyr::n())

  dt_calc_readm <- dt_calc %>%
    dplyr::filter(readmission == TRUE) %>%
    dplyr::group_by(one_month) %>%
    dplyr::summarise(Readm = dplyr::n())


  dt_reamd_disch <- dplyr::left_join(dt_calc_disch, dt_calc_readm) %>%
    tidyr::replace_na(list(Readm = 0))

  # function to plot the 4 hrs emergency performance
  # Plot all days - see Tom's AE APP

  ## Version 0.6.0 of qicharts2 still appears to suffer from
  ## the issue described here: https://github.com/anhoej/qicharts2/issues/21
  ## although it looks like there is a fix. Once that fix is
  ## integrated into a CRAN release, this code should be removed. ***
  options(
    qic.linecol = "#5DA5DA",
    qic.signalcol = "#F15854",
    qic.targetcol = "#059748",
    qic.clshade = TRUE
  )
  ## ***


  # #######################################################
  pct <- qicharts2::qic(Readm,
    n = N,
    x = one_month,
    data = dt_reamd_disch,
    chart = "pp",
    ylab = "percent",
    show.grid = TRUE,
    multiply = 100,
    x.angle = 45
  )

  pct


  # Set the title
  title_stub <- ": Readmissions by "
  start_date_title <- format(as.Date(start_date), format = "%d %B %Y")
  end_date_title <- format(as.Date(end_date), format = "%d %B %Y")
  Days <- " days "
  chart_title <- paste0(hospitalName, title_stub, readmissionBy, Days, start_date_title, " to ", end_date_title)

  pct$data$x <- as.Date(pct$data$x, tz = time_zone)
  cht_data <- add_rule_breaks(pct$data)
  pct <- ggplot2::ggplot(cht_data, ggplot2::aes(x, y, label = x))
  # cutoff <- data.frame(yintercept= 95, cutoff=factor(95))

  # convert arguments to dates and round to nearest quarter
  st.dt <- as.Date(start_date, format = "%Y-%m-%d", tz = time_zone)
  ed.dt <- as.Date(end_date, format = "%Y-%m-%d", tz = time_zone)
  cht_axis_breaks <- seq(st.dt, ed.dt, by = "quarters")
  # ylimlow <- min(min(pct$data$y, na.rm = TRUE),min(pct$data$lcl, na.rm = TRUE))

  readmission_plot <- format_control_chart(pct, r1_col = "orange", r2_col = "steelblue") +
    ggplot2::scale_x_date(
      date_breaks = "1 month", labels = scales::date_format("%Y-%m-%d"),
      breaks = cht_axis_breaks
    ) +
    ggplot2::ggtitle(chart_title) +
    ggplot2::theme(plot.title = ggplot2::element_text(size = 11, face = "bold")) +
    ggplot2::theme(axis.text.x = ggplot2::element_text(angle = 45, size = 10)) +
    ggplot2::labs(
      x = "Month", y = "Percentage of readmissions",
      caption = "*Shewart chart rules apply (see Understanding the Analysis tab for more detail) \nRule 1: Any month outside the control limits \nRule 2: Eight or more consecutive months all above, or all below, the centre line", size = 10
    )

  readmission_plot

  if (returnPlot == TRUE) {
    readmission_plot
  } else {
    readmission_plot$data
  }
}