R/GRscatter.R

In uc-bd2k/GRmetrics_old: Calculate growth-rate inhibition (GR) metrics

#' Scatterplot of a given GR metric
#'
#' Given a SummarizedExperiment object created by \code{\link{GRfit}},
#' this function creates a scatterplot according to the parameters below.
#'
#' @param fitData a SummarizedExperiment object, generated by the GRfit 
#' function.
#' @param metric the GR metric (GR50, GRinf, h_GR, GRmax, GEC50, or GR_AOC) 
#' or traditional metric (IC50, Einf, h, Emax, EC50, or AUC) that will be used 
#' for the scatterplot.
#' @param variable The name of the variable from data (e.g. drug, cell-line,
#' etc.) to select factors from.
#' @param xaxis a vector of values of "variable" of data to be used for the
#' scatterplot's x-axis
#' @param yaxis a vector of values of "variable" of data to be used for the
#' scatterplot's y-axis
#' @param plotly a logical value indicating whether to output a ggplot2 graph
#' or a ggplotly graph
#'
#' @return a ggplot2 or ggplotly scatterplot of the x-axis variable(s) versus
#' the y-axis variable(s) for the given GR metric
#' @author Nicholas Clark
#' @details
#' Given a SummarizedExperiment object created by \code{\link{GRfit}},
#' this function creates a scatterplot of a given GR metric (GR50, GRmax,
#' etc.) or traditional metric (IC50, Emax, etc.) with the "xaxis" values 
#' plotted against the "yaxis" values. The results can be viewed in a static 
#' ggplot image or an interactive plotly graph.
#'
#' The xaxis and yaxis vectors must be of the same length or at least one must
#' be of length one. For each pair of values xaxis[i] and yaxis[i], the
#' function will create a scatterplot (all on one graph) of the specified GR
#' metric. If a vector is of length one, it will be repeated to the length of
#' the other vector.
#' @export
#' @seealso To create the object needed for this function, see
#' \code{\link{GRfit}}. For other visualizations, see \code{\link{GRdrawDRC}}
#' and \code{\link{GRbox}}. For online GR calculator and browser, see
#' \url{http://www.grcalculator.org}.
#' @examples
#' # Load Case A (example 1) input
#' data("inputCaseA")
#' head(inputCaseA)
#' # Run GRfit function with case = "A"
#' output1 = GRfit(inputData = inputCaseA,
#' groupingVariables = c('cell_line','agent', 'perturbation','replicate',
#' 'time'))
#' GRscatter(output1, 'GR50', 'agent', c('drugA','drugD'), 'drugB')
#' GRscatter(output1, 'GR50', 'agent', c('drugA','drugD'), 'drugB',
#' plotly = FALSE)

GRscatter = function(fitData, metric, variable, xaxis, yaxis, plotly = TRUE){
  if(length(xaxis) != length(yaxis)) {
    if(length(xaxis) == 1) {
      xaxis = rep(xaxis, length(yaxis))
    } else if (length(yaxis) == 1) {
      yaxis = rep(yaxis, length(xaxis))
    } else {
      stop('xaxis and yaxis must be of the same length or one must be of
           length 1')
    }
  }
  data = cbind(as.data.frame(SummarizedExperiment::colData(fitData)),
              t(SummarizedExperiment::assay(fitData)))
  if(metric == 'GR50') {
    data$log10GR50 = log10(data$GR50)
    metric = 'log10GR50'
  } else if(metric == 'IC50') {
    data$log10IC50 = log10(data$IC50)
    metric = 'log10IC50'
  } else if(metric == 'h_GR') {
    data$`log2h_GR` = log2(data$h_GR)
    metric = 'log2h_GR'
  } else if(metric == 'h') {
    data$log2h = log2(data$h)
    metric = 'log2h'
  }
  all_data = NULL
  for(i in 1:length(xaxis)) {
    xaxis_data = data[data[[variable]] == xaxis[i],]
    yaxis_data = data[data[[variable]] == yaxis[i],]
    # removing the selected xaxis and yaxis values from the "experiment"
    # column for merging
    temp1 = sub(xaxis[i], '', xaxis_data$experiment)
    temp1 = sub('  ', ' ', temp1)
    xaxis_data$merge = temp1

    temp2 = sub(yaxis[i], '', yaxis_data$experiment)
    temp2 = sub('  ', ' ', temp2)
    yaxis_data$merge = temp2

    merge_data = merge(xaxis_data, yaxis_data, by = 'merge')
    merge_data$cross = as.factor(paste(xaxis[i], 'x', yaxis[i], sep = ' '))
    all_data = rbind(all_data, merge_data)
  }

  x_data = paste0(metric,'.x')
  y_data = paste0(metric,'.y')

  ## Get rid of infinite values
  test_finite_x = which(is.finite( all_data[[x_data]] ))
  test_finite_y = which(is.finite( all_data[[y_data]] ))
  test_finite = intersect(test_finite_x, test_finite_y)
  all_data = all_data[test_finite,]

  padding = 0.05
  scatter_values = c(all_data[[x_data]], all_data[[y_data]])
  ## Get rid of any infinite values
  finite_values = which(is.finite(scatter_values))
  scatter_values = scatter_values[finite_values]
  x_min = min(scatter_values, na.rm = TRUE)
  x_max = max(scatter_values, na.rm = TRUE)
  y_min = min(scatter_values, na.rm = TRUE)
  y_max = max(scatter_values, na.rm = TRUE)
  all_max = max(abs(c(x_max, y_max, x_min, y_min)), na.rm = TRUE)
  all_range = 2*all_max
  all_max = all_max + padding*all_range
  all_min = -all_max

  p = ggplot2::ggplot(data = all_data, ggplot2::aes_string(x = x_data,
                      y = y_data, colour = 'cross', text = 'merge')) +
    ggplot2::geom_point(size=2) +
    ggplot2::geom_abline(slope = 1, intercept = 0, size = .25) +
    ggplot2::coord_fixed() +
    ggplot2::scale_x_continuous(limits = c(all_min, all_max)) +
    ggplot2::scale_y_continuous(limits = c(all_min, all_max)) +
    ggplot2::coord_fixed()

  if(metric == 'log10GR50') {
    p = p + ggplot2::ggtitle("GR50 Scatterplot (log10)") +
      ggplot2::labs(colour = "") + ggplot2::xlab("log10(GR50)") +
      ggplot2::ylab("log10(GR50)")
  } else if(metric == 'log10IC50') {
    p = p + ggplot2::ggtitle("IC50 Scatterplot (log10)") +
      ggplot2::labs(colour = "") + ggplot2::xlab("log10(IC50)") +
      ggplot2::ylab("log10(IC50)")
  } else if(metric == 'log2h_GR') {
    p = p + ggplot2::ggtitle("Hill Slope (h_GR) Scatterplot (log2)") +
      ggplot2::labs(colour = "") + ggplot2::xlab("log2(h_GR)") +
      ggplot2::ylab("log2(h_GR)")
  } else if(metric == 'log2h') {
    p = p + ggplot2::ggtitle("Hill Slope (h) Scatterplot (log2)") +
      ggplot2::labs(colour = "") + ggplot2::xlab("log2(h)") +
      ggplot2::ylab("log2(h)")
  } else {
    p = p + ggplot2::ggtitle(paste(metric, "Scatterplot")) +
      ggplot2::labs(colour = "") + ggplot2::xlab(metric) +
      ggplot2::ylab(metric)
  }
  if(plotly == TRUE) {
    q = plotly::ggplotly(p)
    return(q)
  } else {
    return(p)
  }
}