R/ifcb_psd_plot.R

In iRfcb: Tools for Managing Imaging FlowCytobot (IFCB) Data

utils::globalVariables(c("x", "y"))
#' Generate PSD Plot for a Given Sample
#'
#' This function generates a plot for a given sample from Particle Size Distribution (PSD) data and fits from Imaging FlowCytobot (IFCB).
#' The PSD data and fits can be generated by `ifcb_psd` (Hayashi et al. in prep).
#'
#' @param sample_name The name of the sample to plot in DYYYYMMDDTHHMMSS.
#' @param data A data frame containing the PSD data (data output from `ifcb_psd`), where each row represents a sample and each column represents different particle sizes in micrometers.
#' @param fits A data frame containing the fit parameters for the power curve (fits output from `ifcb_psd`), where each row represents a sample and the columns include the parameters `a`, `k`, and `R2`.
#' @param start_fit The x-value threshold below which data should be excluded from the plot and fit.
#'
#' @return A ggplot object representing the PSD plot for the sample.
#' @export
#'
#' @references Hayashi, K., Walton, J., Lie, A., Smith, J. and Kudela M. Using particle size distribution (PSD) to automate imaging flow cytobot (IFCB) data quality in coastal California, USA. In prep.
#'
#' @seealso \code{\link{ifcb_psd}} \url{https://github.com/kudelalab/PSD}
#'
#' @examples
#' \dontrun{
#' #' # Initialize a python session if not already set up
#' ifcb_py_install()
#'
#' # Analyze PSD
#' psd <- ifcb_psd(feature_folder = 'path/to/features',
#'                 hdr_folder = 'path/to/hdr_data',
#'                 save_data = TRUE,
#'                 output_file = 'psd/svea_2021',
#'                 plot_folder = NULL,
#'                 use_marker = FALSE,
#'                 start_fit = 13,
#'                 r_sqr = 0.5,
#'                 beads = 10 ** 9,
#'                 bubbles = 150,
#'                 incomplete = c(1500, 3),
#'                 missing_cells = 0.7,
#'                 biomass = 1000,
#'                 bloom = 5,
#'                 humidity = NULL)
#'
#' # Plot PSD of the first sample
#' plot <- ifcb_psd_plot(sample_name = "D20230316T101514",
#'                       data = psd$data,
#'                       fits = psd$fits,
#'                       start_fit = 10)
#'
#' # Inspect plot
#' print(plot)
#' }
ifcb_psd_plot <- function(sample_name, data, fits, start_fit) {
  # Extract the sample data
  sample_data <- data %>% filter(sample == sample_name)
  x_values <- as.numeric(gsub("[^0-9]+", "", colnames(sample_data)[4:ncol(sample_data)]))
  y_values <- as.numeric(sample_data[1, 4:ncol(sample_data)])

  # Create a data frame for plotting and fitting
  plot_data <- data.frame(x = x_values, y = y_values)

  # Filter out values below start_fit
  plot_data <- plot_data %>% filter(x >= start_fit)

  # Extract the fit parameters
  fit_params <- fits %>% filter(sample == sample_name)

  if (nrow(fit_params) == 0) {
    stop("No fit parameters found for the specified sample.")
  }

  a <- fit_params$a
  k <- fit_params$k
  R2 <- fit_params$R.2

  # Create the equation text
  if (!is.na(R2) && is.finite(R2)) {
    equation_text <- paste0("y = ", format(a, scientific = TRUE, digits = 3), " * x^", format(k, digits = 3), "\nR\u00B2 = ", format(R2, digits = 3))
  } else {
    equation_text <- "No R\u00B2 value available."
  }

  # Plot the data
  p <- ggplot(plot_data, aes(x = x, y = y)) +
    geom_line() +
    labs(title = paste("Sample:", sample_name),
         x = "ESD (\u00B5m)",
         y = "N'(D) [c/L\u207B]") +
    theme_minimal() +
    theme(
      panel.grid.major = element_blank(),
      panel.grid.minor = element_blank(),
      panel.border = element_rect(color = "black", fill = NA),
      panel.background = element_rect(fill = "white", color = NA),
      plot.background = element_rect(fill = "white", color = NA)
    )

  # Add the power curve fit if R2 is not -Inf
  if (!is.na(R2) && is.finite(R2)) {
    p <- p + stat_function(fun = function(x) { a * x^k }, color = "blue") +
      annotate("text", x = max(plot_data$x) * 0.5, y = max(plot_data$y) * 0.9, label = equation_text, hjust = 0, vjust = 1, size = 5, color = "black")
  }

  p
}