R/plot_beta_fcm.R

In rprops/Phenoflow_package: Advanced analysis of microbial flow cytometry data

#' Plot function for beta diversity analysis of FCM data
#'
#' This function visualizes the beta diversity analysis generated by beta_div_fcm()
#' @param x flowbasis object generated by flowBasis()
#' @param color A vector of factors indicating the groups to be colored in the plot
#' @param shape A vector of factors indicating the shape of the groups in the plot
#' @param labels A vector of length 2, giving the labels for color and shape to be 
#' used in the legend. Has to be in the order c('color','shape'). 
#' @param legend.pres Surpresses legend when there is no definition of labels.
#' @keywords diversity, fcm, beta
#' @examples
#' ## Short example
#' 
#' # Load precomputed fingerprint object
#' data(CoolingTower)
#' 
#' # Calculate diversity values
#' beta <- beta_div_fcm(CoolingTower,ord.type="PCoA")
#' plot_beta_fcm(beta)
#' 
#' ## Full data processing example
#' 
#' # Load raw data (imported using flowCore)
#' data(flowData)
#' # Asinh transform and select parameters of interest (cells were stained with Sybr Green I).
#' flowData_transformed <- flowCore::transform(flowData,`FL1-H`=asinh(`FL1-H`),
#'        `SSC-H`=asinh(`SSC-H`), 
#'        `FL3-H`=asinh(`FL3-H`), 
#'        `FSC-H`=asinh(`FSC-H`))
#' param=c('FL1-H', 'FL3-H','SSC-H','FSC-H')
#' flowData_transformed = flowData_transformed[,param]
#' 
#' # Create a PolygonGate for denoising the dataset
#' # Define coordinates for gate in sqrcut1 in format: c(x,x,x,x,y,y,y,y)
#' sqrcut1 <- matrix(c(8.75,8.75,14,14,3,7.5,14,3),ncol=2, nrow=4)
#' colnames(sqrcut1) <- c('FL1-H','FL3-H')
#' polyGate1 <- flowCore::polygonGate(.gate=sqrcut1, filterId = 'Total Cells')
#' 
#' # Gating quality check
#' flowViz::xyplot(`FL3-H` ~ `FL1-H`, data=flowData_transformed[1], filter=polyGate1,
#'          scales=list(y=list(limits=c(0,14)),
#'          x=list(limits=c(6,16))),
#'          axis = lattice::axis.default, nbin=125, 
#'          par.strip.text=list(col='white', font=2, cex=2), smooth=FALSE)
#'  
#'  # Isolate only the cellular information based on the polyGate1
#'  flowData_transformed <- flowCore::Subset(flowData_transformed, polyGate1)
#'  
#'  # Normalize parameter values to [0,1] interval based on max. value across parameters
#'  summary <- flowCore::fsApply(x=flowData_transformed,FUN=function(x) apply(x,2,max),use.exprs=TRUE)
#'  max = max(summary[,1])
#'  mytrans <- function(x) x/max
#'  flowData_transformed <- flowCore::transform(flowData_transformed,`FL1-H`=mytrans(`FL1-H`),
#'          `FL3-H`=mytrans(`FL3-H`), 
#'          `SSC-H`=mytrans(`SSC-H`),
#'          `FSC-H`=mytrans(`FSC-H`))
#'  
#'  # Calculate fingerprint
#'  fbasis <- flowFDA::flowBasis(flowData_transformed, param, nbin=128, 
#'          bw=0.01, normalize=function(x) x)
#'  
#'  # Calculate diversity
#'  beta <- beta_div_fcm(fbasis,ord.type="PCoA")
#'  plot_beta_fcm(beta)

#' @export

plot_beta_fcm <- function(x, color = NA, shape = NA, labels = c("Factor 1", "Factor 2"), legend.pres = NULL) {
  legend.ops <- NULL
  if (sum(is.na(color)) > 0) 
    color = rep("f1", nrow(x$points))
  if (sum(is.na(shape)) > 0) {
    shape = rep("f2", nrow(x$points))
    legend.ops <- FALSE
  }
  var.pcoa <- vegan::eigenvals(x)/sum(vegan::eigenvals(x))
  PcoA <- as.data.frame(x$points)
  names(PcoA)[1:2] <- c("Axis1", "Axis2")
  PcoA <- cbind(PcoA, color, shape)
  ggplot2::ggplot(PcoA, ggplot2::aes(x = Axis1, y = Axis2, color = color, shape = shape)) + 
    ggplot2::geom_point(alpha = 0.7, size = 4) + 
    ggplot2::geom_point(colour = "grey90", size = 1.5) +
    ggplot2::scale_color_manual(values = c("#a65628", "red", "#ffae19",
                                  "#4daf4a", "#1919ff", "darkorchid3", "magenta")) +
    ggplot2::labs(x = paste0("Axis1 (",round(100*var.pcoa[1],1),"%)"),
                  y = paste0("Axis2 (",round(100*var.pcoa[2],1),"%)")) +
    ggplot2::ggtitle("Ordination of phenotypic fingerprints") +
    ggplot2::labs(color = labels[1],shape = labels[2]) +
    ggplot2::guides(color = legend.pres, shape = legend.ops)
}