R/GRdrawDRC.R

In GRmetrics: Calculate growth-rate inhibition (GR) metrics

#' Dose-Response Curves
#'
#' Given a SummarizedExperiment object created by \code{\link{GRfit}},
#' this function plots either the growth-rate inhibition (GR) dose response 
#' curves or the traditional dose response curves for a given set of data.
#'
#' @param fitData an element of class SummarizedExperiment, generated by the 
#' GRfit function.
#' @param metric either "GR" for GR dose response curves or "rel_cell" for 
#' traditional dose response curves based on relative cell count.
#' @param experiments the names of the experiments to plot (or "all")
#' @param min the minimum concentration to plot (for curves)
#' @param max the maximum concentration to plot (for curves)
#' @param points a logical value indicating whether points (individual GR
#' values) will be plotted
#' @param curves a logical value indicating whether sigmoidal dose-response
#' curves will be plotted
#' @param plotly a logical value indicating whether to output a ggplot2 graph
#' or a ggplotly graph
#'
#' @return ggplot2 or ggplotly graphs of Growth-rate inhibition dose-response
#' curves
#' @author Nicholas Clark
#' @details
#' Given a SummarizedExperiment object created by \code{\link{GRfit}},
#' this function plots these GR values (versus concentration) and/or the
#' sigmoidal curves fitted to the sets of points. The results can be viewed
#' in a static ggplot image or an interactive plotly graph.
#'
#' The "min" and "max" parameters control the concentration values for which
#' the curves are plotted. They are automatically set to the minimum and
#' maximum concentrations of the data, but can be set by the user as well.
#' "min" and "max" take raw values (not log transformed) for concentration.
#'
#' By default, curves and points are plotted for all experiments. To specify
#' a smaller set of experiments, use the "experiments" parameter. To see the
#' names of individual experiments for a GRfit object \code{fit_example}, see
#' \code{colData(fit_example)}. See the examples below.
#' @seealso To create the object needed for this function, see
#' \code{\link{GRfit}}. For other visualizations, see \code{\link{GRbox}} and
#' \code{\link{GRscatter}}. For online GR calculator and browser, see
#' \url{http://www.grcalculator.org}.
#' @examples
#' # Load Case A (example 1) input
#' data("inputCaseA")
#' # Run GRfit function with case = "A"
#' drc_output = GRfit(inputCaseA,
#' groupingVariables = c('cell_line','agent'))
#' GRdrawDRC(drc_output, experiments = c('BT20 drugA', 'MCF10A drugA',
#' 'MCF7 drugA'), min = 10^(-4), max = 10^2)
#' GRdrawDRC(drc_output, plotly = FALSE)
#' @export

GRdrawDRC <- function(fitData, metric = "GR", experiments = "all",
                      min = "auto", max = "auto",
                      points = TRUE, curves = TRUE, plotly = TRUE) {
  if(points == FALSE & curves == FALSE) {
    stop('You must show either points or curves or both')
  }
  if(metric == "IC") {
    message('For the traditional dose-response curve based on relative cell counts, please use metric = "rel_cell" instead of metric = "IC". This notation has been changed as of Version 1.3.2.')
  }
  # declaring values NULL to avoid note on package check
  log10_concentration = NULL
  experiment = NULL
  data = S4Vectors::metadata(fitData)[[1]]
  parameterTable = cbind(as.data.frame(SummarizedExperiment::colData(fitData)),
                        t(SummarizedExperiment::assay(fitData)))
  groupingVariables = S4Vectors::metadata(fitData)[[2]]
  data$log10_concentration = log10(data$concentration)
  tmp<-data[,groupingVariables, drop = FALSE]
  experimentNew = (apply(tmp,1, function(x) (paste(x,collapse=" "))))
  if(length(groupingVariables) > 0) {
    data$experiment = as.factor(experimentNew)
  } else {
    data$experiment = as.factor("All Data")
  }
  if(!identical(experiments, "all")) {
    parameterTable = parameterTable[parameterTable$experiment %in%
                                      experiments, ]
    data = data[data$experiment %in% experiments, ]
  }
  exps = unique(parameterTable$experiment)
  if(min == "auto") {
    min_conc = min(data$concentration, na.rm = TRUE)
  } else {
    min_conc = min
  }
  if(max == "auto") {
    max_conc = max(data$concentration, na.rm = TRUE)
  } else {
    max_conc = max
  }
  len = (log10(max_conc) - log10(min_conc))*100
  Concentration = 10^(seq(log10(min_conc) - 1, log10(max_conc) + 1,
                          length.out = len))
  curve_data_all = NULL
  for(exp in exps) {
    row = which(parameterTable$experiment == exp)
    if(metric == "GR") {
      GEC50 = parameterTable$GEC50[row]
      GRinf = parameterTable$GRinf[row]
      h_GR = parameterTable$h_GR[row]
      logistic_3u = function(c){GRinf + (1 - GRinf)/(1 + (c/GEC50)^h_GR)}
    } else if (metric %in% c("rel_cell", "IC")) {
      EC50 = parameterTable$EC50[row]
      Einf = parameterTable$Einf[row]
      h = parameterTable$h[row]
      logistic_3u = function(c){Einf + (1 - Einf)/(1 + (c/EC50)^h)}
    }
    curve_data = as.matrix(Concentration)
    colnames(curve_data) = "Concentration"
    if(metric == "GR") {
      if(parameterTable$fit_GR[row] == "sigmoid") {
        GR = apply(curve_data, 1, logistic_3u)
      } else {
        GR = parameterTable$flat_fit_GR[row]
      }
      curve_data = cbind(curve_data, GR)
    } else if(metric %in% c("rel_cell", "IC")) {
      if(parameterTable$fit_rel_cell[row] == "sigmoid") {
        rel_cell_count = apply(curve_data, 1, logistic_3u)
      } else {
        rel_cell_count = parameterTable$flat_fit_rel_cell[row]
      }
      curve_data = cbind(curve_data, rel_cell_count)
    }
    curve_data = as.data.frame(curve_data)
    curve_data$experiment = exp
    if(is.null(curve_data_all)){
      curve_data_all = curve_data
    } else {
      curve_data_all = rbind(curve_data_all, curve_data)
    }
  }
  curve_data_all = merge(curve_data_all, parameterTable[, c(groupingVariables, "experiment")])
  curve_data_all$experiment = as.factor(curve_data_all$experiment)

  if(metric == "GR") {
    if(points == TRUE & curves == FALSE) {
      p = ggplot2::ggplot(data = data, ggplot2::aes(x = log10_concentration,
                          y = GRvalue, colour = experiment)) + ggplot2::geom_point()
    } else if(points == FALSE & curves == TRUE) {
      p = ggplot2::ggplot(data = curve_data_all,
                          ggplot2::aes(x = log10(Concentration),
                          y = GR, colour = experiment)) + ggplot2::geom_line()
    } else if(points == TRUE & curves == TRUE) {
      p = ggplot2::ggplot() + ggplot2::geom_line(data = curve_data_all,
        ggplot2::aes(x = log10(Concentration), y = GR,colour = experiment)) +
        ggplot2::geom_point(data = data, ggplot2::aes(x = log10_concentration,
                                              y = GRvalue, colour = experiment))
    }
    p = p + ggplot2::coord_cartesian(xlim = c(log10(min_conc)-0.1,
                                              log10(max_conc)+0.1),
                                     ylim = c(-1, 1.5), expand = TRUE) +
      ggplot2::ggtitle("Concentration vs. GR values") +
      ggplot2::xlab('Concentration (log10 scale)') +
      ggplot2::ylab('GR value') + ggplot2::labs(colour = "") +
      ggplot2::geom_hline(yintercept = 1, size = .25) +
      ggplot2::geom_hline(yintercept = 0.5, size = .25) +
      ggplot2::geom_hline(yintercept = 0, size = .25) +
      ggplot2::geom_hline(yintercept = -1, size = .25)
  } else if(metric %in% c("rel_cell", "IC")) {
    if(points == TRUE & curves == FALSE) {
      p = ggplot2::ggplot(data = data, ggplot2::aes(x = log10_concentration,
          y = rel_cell_count, colour = experiment)) + ggplot2::geom_point()
    } else if(points == FALSE & curves == TRUE) {
      p = ggplot2::ggplot(data = curve_data_all,
                          ggplot2::aes(x = log10(Concentration),
          y = rel_cell_count, colour = experiment)) + ggplot2::geom_line()
    } else if(points == TRUE & curves == TRUE) {
      p = ggplot2::ggplot() + ggplot2::geom_line(data = curve_data_all,
        ggplot2::aes(x = log10(Concentration), y = rel_cell_count,
                     colour = experiment)) +
        ggplot2::geom_point(data = data, ggplot2::aes(x = log10_concentration,
                            y = rel_cell_count, colour = experiment))
    }
    p = p + ggplot2::coord_cartesian(xlim = c(log10(min_conc)-0.1,
                                              log10(max_conc)+0.1),
                                     ylim = c(0, 1.5), expand = TRUE) +
      ggplot2::ggtitle("Concentration vs. Relative cell count") +
      ggplot2::xlab('Concentration (log10 scale)') +
      ggplot2::ylab('Relative cell count') + ggplot2::labs(colour = "") +
      ggplot2::geom_hline(yintercept = 1, size = .25) +
      ggplot2::geom_hline(yintercept = 0.5, size = .25) +
      ggplot2::geom_hline(yintercept = 0, size = .25)
  }
  if(plotly == TRUE) {
    q = plotly::ggplotly(p)
    return(q)
  } else {
    return(p)
  }
}