R/estimateEdgeEffect.R

In lujunyan1118/DrugScreenExplorer_dev: A pipeline for processing drug sensitivity screen data

#' Estimate edge effect intensity
#'
#' This function estimates the intensity of edge effect on each screen plate by comparing the median signal intensities on the edge to those on the inner plate.
#' The results will be shown as a bar plot with plate IDs as x-axis and edge effect intensities as y-axis.
#'
#' @param screenData a data frame containing the screening results generated by \code{readScreen()} function
#' @param nLayer an integer value (larger than 1), which indicates up to how many layers on the edge of the screen are considered as edges. The default value is 1.
#' @param onlyNeg a logical value, whether only wells annotated as negative controls are used for edge effect estimation.
#' If \code{FALSE}, all the measurements on each plate are used for edge effect estimation. The default value is \code{TRUE}.
#' @param pdfName either NULL or a character string specifying the output file name for the plot in pdf format. If NULL, no plot will be created.
#' @param identifier a character string specifying the column in the input data frame that should be used as identifiers (x-axis) in the plot.
#' Default value is "fileName". If one identifier matches several plates, error bars indicating standard deviations will be shown in the plot.
#' @export
#' @import ggplot2 dplyr
#' @return Depends on the number \emph{N} specified with the \code{nLayer} argument, this function will add \emph{N} columns with the names edgeRatio1, edgeRatio2, ... edgeRation\emph{N}. Each column contains the estimated edge effect intensity when considering the \emph{N} outer layers on the plate as edge layers.
#' If a file name is specified for the \code{pdfName} argument, a pdf file that contains a bar plot visualizing the edge effect intensity will be created in the working directly.
#' @examples
#' # load processed data
#' data('screenData_normalized')
#'
#' # estimate the intensity of edge effect.
#' # The wells from the first outside layer on the plate will be defined as the edge wells.
#' estimateEdgeEffect(screenData_normalized, nLayer =1)
#'
#' # Please see the vignette for more information.

estimateEdgeEffect <- function(screenData, nLayer = 1, onlyNeg = TRUE, identifier = "fileName",
                               pdfName = NULL) {

  for (n in seq(nLayer)) {
    varName <- paste0("edgeRatio", n)

    edgeEffects <- lapply(unique(screenData$fileName), function(eachFile) {
      screenSub <- dplyr::filter(screenData, fileName %in% eachFile)
      outTable <- tibble(fileName = eachFile, val = estimateOnePlate(screenData,
                                                                     n, onlyNeg))
      colnames(outTable) <- c("fileName", varName)
      outTable

    }) %>%
      bind_rows()

    screenData <- left_join(screenData, edgeEffects, by = "fileName")
  }

  if (!is.null(pdfName)) {
    # prepare table for plots
    mapCol <- structure(seq(nLayer), names = paste0("edgeRatio", seq(nLayer)))
    plotTab <- rename(screenData, idn = !!identifier) %>%
      dplyr::select(idn, names(mapCol)) %>%
      gather(key = "layer", value = "edgeRatio", -idn) %>%
      dplyr::mutate(layer = factor(mapCol[layer])) %>%
      group_by(idn, layer) %>%
      summarise(meanFac = mean(edgeRatio), sdFac = sd(edgeRatio)) %>%
      arrange(desc(meanFac)) %>%
      ungroup() %>%
      dplyr::mutate(idn = factor(idn, levels = unique(idn)))

    # create plot
    p <- ggplot(plotTab, aes(x = idn, y = meanFac, fill = layer)) +
      geom_bar(position = position_dodge(),stat = "identity") +
      coord_flip() + theme_bw() + xlab(identifier) + ylab("edge effect ratio")

    # generate pdf
    pdfHeight <- length(unique(plotTab$idn)) * 0.8
    pdf(file = pdfName, width = 10, height = pdfHeight)
    plot(p)
    dev.off()
  }

  screenData
}