Microbiome
McConville Lab
BMDM Subpopulation Analysis
Condition: Unstimulated Plate A

workflowr

• Summary
• Checks
• Past versions

Last updated: 2025-03-24

Checks: 7 0

Knit directory: McConville_Lab_Microbiome_PMC255A/

This reproducible R Markdown analysis was created with workflowr (version 1.7.1). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.

---

R Markdown file: up-to-date

Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.

Environment: empty

Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.

Seed: set.seed(20250117)

The command set.seed(20250117) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.

Session information: recorded

Great job! Recording the operating system, R version, and package versions is critical for reproducibility.

Cache: none

Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.

File paths: relative

Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.

Repository version: 3a65c0b

Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.

The results in this page were generated with repository version 3a65c0b. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.

Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:

Ignored files:
    Ignored:    .Rhistory
    Ignored:    .Rproj.user/
    Ignored:    data/raw_data/

Untracked files:
    Untracked:  analysis/PMC255A_ScreenQC_v4.Rmd
    Untracked:  analysis/prelim/

Unstaged changes:
    Modified:   analysis/PMC255A_ScreenComparison.Rmd
    Deleted:    analysis/PMC255A_ScreenQC_v3.Rmd
    Modified:   analysis/_site.yml

Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.

---

These are the previous versions of the repository in which changes were made to the R Markdown (analysis/PMC255A_Unstim_PlateA_DT.Rmd) and HTML (docs/PMC255A_Unstim_PlateA_DT.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.

File	Version	Author	Date	Message
Rmd	4df8175	annrann.wong	2025-02-13	Update Viability Bins
html	4df8175	annrann.wong	2025-02-13	Update Viability Bins
Rmd	fc42795	annrann.wong	2025-02-12	Update MeanSpotCount readout
html	fc42795	annrann.wong	2025-02-12	Update MeanSpotCount readout
Rmd	25b8bb6	annrann.wong	2025-02-12	Update Analysed Data Table
html	25b8bb6	annrann.wong	2025-02-12	Update Analysed Data Table

---

Aim

To examine the polarisation potential of BDBM cells from M0 to M1 or M2 in the presence of compounds and different concentrations of polarisation activator.

Researcher

Ada Koo, McConville Lab, Bedoui Lab

Readout

Image analysis of BMDM cells stained with DAPI, ARG, INOS, and BODIPY (CX7 Pro LZR, 9 fields0 @ 20x)

Image analysis software

CellProfiler 4.1.3

Required R packages

data.table, DT, reshape2, tidyverse, viridis, patchwork

# load the required packages
library(data.table)
library(readxl)
library(DT)
library(reshape2)
library(viridis)
library(patchwork)
library(tidyverse)
library(magick)

# fix ggplotly axis labels
layout_ggplotly <- function(gg, x = -0.02, y = -0.08){
  # The 1 and 2 goes into the list that contains the options for the x and y axis labels respectively
  gg[['x']][['layout']][['annotations']][[1]][['y']] <- x
  gg[['x']][['layout']][['annotations']][[2]][['x']] <- y
  gg
}

# set the file prefix
prefix <- "KW_PMC255A"

# turn off scientific notation
options(scipen=999)

---

Screen details

Screen dates (yyyy-mm-dd)

2024-12-04

Cell lines

BMDM M0

Conditions

unstimulated (Unstim_M0)
5 ng/mL IL4 + 5 ng/mL IL13 (Lo_M2)
10 ng/mL IL4 + 10 ng/mL IL13 (Hi_M2)
3 ng/mL IFNg + 3 ng/mL LPS (Lo_M1)
25 ng/mL IFNg + 25 ng/mL LPS (Hi_M1)
25 ng/mL IFNg (Potential_M1)

Screen outline

Library: Custom MCE library
Library plates: 4
Cell plates per library plate: 24

Controls

Per well data is normalised to Stimulated DMSO on a per plate basis.

---

QC considerations

None.

---

Plate layout

get_image_dimensions <- function(path) {
  img_info <- magick::image_info(magick::image_read(path))
  setNames(c(img_info$width, img_info$height), c("Width", "Height"))
}

# Example usage
path <- file.path("data", paste0(prefix, "_platemap.png"))
img_dims <- get_image_dimensions(path)

# Generate the HTML image with scaled dimensions
htmltools::img(
  src = knitr::image_uri(path),
  style = sprintf("height:%.0fpx;width:%.0fpx;", img_dims["Height"] / 2, img_dims["Width"] / 2)
)

---

Data

The raw and normalised (fold change to DMSO) data for all wells can be downloaded below.

data_norm <- read_csv(paste0("output/", prefix, "_data_norm.csv")) %>% 
  filter(grepl("Unstim_M0", Barcode))
data_norm_melt <- read_csv(paste0("output/", prefix, "_data_norm_melt.csv")) %>% 
  filter(grepl("Unstim_M0", Barcode))

# print cells in table
condition_list <- unique(data_norm$Condition_Name)

  for (condition in condition_list){
    cat("\n##", condition, "{.tabset .tabset-fade .tabset-pills}\n")  
  
  data_norm_filt <- data_norm %>%
  filter(Condition_Name == condition) %>% 
  select(-Condition_Name) %>% 
  mutate_if(is.numeric, round, 2)

# display the data in a table
 print(htmltools::tagList(datatable(data_norm_filt, extensions = c('Buttons', 'Scroller'), rownames = FALSE, 
  options = list(
  dom = 'Blfrtip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
  pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "_Norm_Data_", condition)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "_Norm_Data_", condition), 
                      title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))
}

Unstimulated

---

Data processing

Robust Z-Scoring

The normalised values have been Robust Z-Scored to the median of all of the test compounds in order to assess the relative strength of each compound. This method can be used to select the strongest and most robust hits from a large dataset. Z-Scores of <= -2 or >= 2 are considered to be significant, but this value must be viewed in combination with the fold change. See the Robust Z-Scoring section of the Methods page for a more information.

# calculate the dmso med and mad
data_zscored_melt <- data_norm_melt %>%
  filter(Description == "Sample") %>% # & !grepl("Plate4", Barcode)
  group_by(Feature) %>%
  mutate(Robust_ZScore = (NormToStimDMSO - median(NormToStimDMSO, na.rm = TRUE))/mad(NormToStimDMSO, na.rm = TRUE)) %>%
  mutate_if(is.numeric, round, 2) %>% 
  select(-Feature, -Raw, -NormToUnstimDMSO, -NormToStimDMSO, -Robust_ZScore, -Phenotype, -Stimulation, Feature, Raw, NormToStimDMSO, Robust_ZScore)

# "unmelt" the dataframe
data_zscored <- data_zscored_melt %>% 
  pivot_wider(id_cols = colnames(data_zscored_melt[c(1:12)]), names_from = Feature, values_from = c(Raw, NormToStimDMSO, Robust_ZScore),
              names_glue = "{.value}_{Feature}")

Viability binning

Compounds were sorted into the following viability bins based on the normalised DMSO (stimulated):

• High viability (> 1.15)
• Normal viability (0.8 <= 1.15)
• Moderate viability (0.5 <= 0.8)
• Low viability (<= 0.5)
• Very Low viability (<= 0.3)

# add viability bins to data
analysed_data_melt <- data_zscored_melt %>% 
    filter(!grepl("Positive_Perc_BODIPY|_PerBODIPYPosCell", Feature)) %>% 
  mutate(Viability_Bin = 
           case_when(grepl("DAPI", Feature) & NormToStimDMSO > 1.15 ~ "High",
                     grepl("DAPI", Feature) & NormToStimDMSO > 0.8 & NormToStimDMSO <= 1.15 ~ "Normal",
                     grepl("DAPI", Feature) & NormToStimDMSO > 0.5 & NormToStimDMSO <= 0.8 ~ "Moderate",
                     grepl("DAPI", Feature) & NormToStimDMSO > 0.3 & NormToStimDMSO <= 0.5 ~ "Low",
                     grepl("DAPI", Feature) & NormToStimDMSO <= 0.3 ~ "VeryLow")) 

analysed_data <- analysed_data_melt %>% 
  pivot_wider(id_cols = c(Barcode, QCL_Sample_Number, Compound_Name, Concentration, Units, WellID, Target, Pathway), names_from = c(Feature), values_from = matches("^Raw|^NormToStimDMSO|^Robust_ZScore|Viability_Bin")) %>% 
  select(-contains(c("Viability_Bin_Positive_Perc", "Viability_Bin_MeanSpotCount"))) %>% 
  mutate_if(is.numeric, round, 2)

---

Analysed Data

The analysed data for all test compounds can be downloaded below.

    print(htmltools::tagList(datatable(analysed_data, extensions = c("Buttons", "Scroller", "FixedColumns", "ColReorder"), rownames = FALSE,
      options = list(
        dom = "Bflrtip",
        columnDefs = list(list(className = "dt-center", targets = "_all")),
        lengthMenu = list(c(5, -1), c("5", "All")),
        pageLength = 5,
        buttons = list("colvis",
          list(extend = "csv", filename = paste0(prefix, "_AnalysedData_", condition)),
          list(extend = "excel", filename = paste0(prefix, "_AnalysedData_", condition), title = NULL)),
        scrollX = TRUE, fixedColumns = list(leftColumns = 2)))))

  cat("\n")

ARG Perc Positive

Waterfall plots

Robust Z-Scores for all of the test compounds are displayed below. The plots is interactive. Hover over the points to view the compound name and other information. The number of hits (unique compound IDs) and corresponding fold change (FC) at Z-Score cut-offs of 2, 3 and 4 are printed below the plots.

# prepare the data
feat_list <- unique(data_zscored_melt$Feature)[2] #"Positive_Perc_ARG"

data_zscored_plot <- data_zscored_melt %>%
  filter(Feature == feat_list) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  unite(Treatment, c(Compound_Name, Conc), remove = FALSE) %>%
  mutate(Status = ifelse(Robust_ZScore >= 2, "Significant Increase",
                  ifelse(Robust_ZScore <= -2, "Significant Decrease", "No Change")))

# set order of status for plotting
stat_order <- c("Significant Increase", "Significant Decrease", "No Change")
data_zscored_plot$Status <- factor(data_zscored_plot$Status, levels = stat_order)

# reorder
data_zscored_plot_filt <- data_zscored_plot %>%
  group_by(Feature) %>% 
  arrange(Robust_ZScore) %>%
  mutate(Index = 1:n())

# print a plot and dataframe for each feature
# for (condition in cond_list){
  
# create the condition heading
# cat("\n")
# cat("\n### ", condition, "\n")
# cat("\n")

  for (feature in feat_list) {
  
  # create the feature heading
  cat("\n")
  cat("\n#### ", feature, "\n")
  cat("\n")

# create dot plot
snake_plot <- data_zscored_plot_filt %>%
  filter(Feature == feature) %>% 
  ggplot(., aes(x = Index, 
                y = Robust_ZScore,
                label = Compound_Name)) +
  geom_point(aes(colour = Status,
                 text = paste('Compound_Name:', Compound_Name,
                              '<br>Raw:', Raw,
                              '<br>NormToStimDMSO:', NormToStimDMSO,
                              '<br>Robust Z-Score:', Robust_ZScore,
                              '<br>PlateID:', PlateID)), 
                 alpha = 0.8, size = 1.6) +
  geom_hline(aes(yintercept = -2), linetype = "dotted", size = 0.3) +
  geom_hline(aes(yintercept = 2), linetype = "dotted", size = 0.3) +
  labs(x = "Compound",
       y = "Z-Score To Sample") +
  #scale_y_continuous(breaks = seq(-6, 6, by = 1), limits=c(-6,6)) +
  theme_bw() +
  theme(#axis.title.y = element_text(margin = margin(t = 0, r = 30, b = 0, l = 30)), 
        axis.title.y = element_text(margin = margin(t = 0, r = 10, b = 0, l = 10)),
        axis.title.x = element_text(size = 1, colour = "white"),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_text(colour = "black"),
        plot.title = element_text(margin = margin(t = 10, r = 0, b = 5, l = 0), size = 14), 
        plot.margin = margin(t = 10, r = 30, b = 5, l = 30),
        legend.title = element_blank()) +
  scale_colour_manual(values = c("#bc3754", "#2a788e", "grey61"), drop = FALSE, name = "")

print(htmltools::tagList(plotly::ggplotly(snake_plot, tooltip = "text")))

cat("\n")
cat("\n")

data_zscored_filt <- data_zscored_melt %>%
  filter(Feature == feature)

hit_count_summary_2 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_3 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 3, ">= 3",
                      ifelse(Robust_ZScore <= -3, "<= -3",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_4 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 4, ">= 4",
                      ifelse(Robust_ZScore <= -4, "<= -4",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary <- rbind(hit_count_summary_2, hit_count_summary_3, hit_count_summary_4) %>% 
  pivot_longer(c(Equivalent_FC_max, Equivalent_FC_min), 
               names_to = "Direction", values_to = "Equivalent_FC") %>% 
  filter(c(grepl('<', ZScore_Bin) & Direction == "Equivalent_FC_max") |
         c(grepl('>', ZScore_Bin) & Direction == "Equivalent_FC_min")) %>% 
  mutate(ZScore_Bin = factor(ZScore_Bin, levels = c(">= 2", ">= 3", ">= 4", "<= -2", "<= -3", "<= -4"))) %>%
  complete(ZScore_Bin, fill = list(Hit_Count = 0)) %>% 
  select(ZScore_Bin, Equivalent_FC, Hit_Count) %>%
  arrange(ZScore_Bin) %>% 
  mutate_if(is.numeric, round, 2) 
  
# display the table
print(htmltools::tagList(datatable(hit_count_summary, extensions = c('Scroller', "FixedColumns"), rownames = FALSE, 
  options = list(
  dom = 'rt',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
          pageLength = -1,
  searching = FALSE,
  scrollX = TRUE))))

  cat("\n")
  } # end the loop

Positive_Perc_ARG

# }

# Create vectors for WellID, row, and column
WellID <- character(384)
row <- character(384)
column <- integer(384)
# list_of_celllines <- unique(sort(data_zscored$Cell_Line))

# Generate WellID, row, and column values
for (i in 1:384) {
  row[i] <- LETTERS[ceiling(i/24)] # Get the row letter
  column[i] <- i %% 24              # Get the column number
  if (column[i] == 0) column[i] <- 24 # Adjust column number for multiples of 24
  WellID[i] <- paste0(row[i], sprintf("%02d", column[i])) # Combine row and column to form WellID
}

# Create a data frame
row_labels <- data.frame(WellID = WellID, row = row, column = column)

# define factor levels
counts <- c(0:9, "10+")

# function to define colour palette and matching values
scale_fill_hitcount <- function(...){
    ggplot2:::manual_scale(
        'fill',
        values = setNames(
          c("#440154", "#482878", "#3e4989", "#31688e", "#26828e",
            "#1f9e89", "#35b779", "#6ece58", "#b5de2b", "#fde725", "#f98e09"),
          counts),
        ...
    )
}

# count the plates
  plate_count <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    pull(VCFG_Plate_ID) %>%
    unique() %>%
    length()

  # count the hits per well/run
  hit_count_summary_2_well <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                        ifelse(Robust_ZScore <= -2, "<= -2",
                        "NonSig"))) %>%
    filter(!ZScore_Bin == "NonSig") %>%
    group_by(WellID, ZScore_Bin) %>%
    summarise(Hit_Count = n())

  # zscore increase
  data_plot1 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == ">= 2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = ">= 2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # zscore decrease
  data_plot2 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == "<= -2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = "<= -2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # make plot
  plot1 <- data_plot1 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    ggtitle(paste("Across", plate_count, "plate\n")) +
    labs(subtitle = " Robust_ZScore >= 2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  # make plot
  plot2 <- data_plot2 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    labs(subtitle = " Robust_ZScore <= -2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

cols <- c("#23888e", "royalblue1", "darkblue")

# count the hits per well/run
hit_count_summary_2_plate <- data_zscored_melt %>%
  filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig")),
         VCFG_Plate_ID = formatC(VCFG_Plate_ID, width = 2, format = "d", flag = "0")) %>%
  unite(VCFG_Plate_ID, c(Barcode), remove = FALSE) %>%
  group_by(VCFG_Plate_ID, Barcode, Conc, ZScore_Bin) %>%
  summarise(Hit_Count = n()) 

# create a list of plates/conc
all_plates <- distinct(hit_count_summary_2_plate[1:3])


# zscore increase
data_plot1 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == ">= 2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = ">= 2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# zscore decrease
data_plot2 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == "<= -2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = "<= -2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# make plot
plot1 <- data_plot1 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore >= 2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  # ylim(0, 15) +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_blank(),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.4, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 0, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

plot2 <- data_plot2 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore <= -2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_text(size = 11, margin = margin(t = 10, r = 0, b = 0, l = 0)),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.5, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 10, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

# } # end the loop

Hit selection ARG (Z-score)

Compounds with a negative Z score are ARG % positive increasing hits while compounds with a positive Z score are ARG % positive decreasing hits.

# include DMSO percentage positive wells
  dmsopercpositive <- data_norm %>% 
    filter(Compound_Name == "DMSO" & !grepl("01", WellID)) %>%
    group_by(Compound_Name) %>% 
    summarise(arg_pos_perc = median(Raw_Positive_Perc_ARG),
              inos_pos_perc = median(Raw_Positive_Perc_INOS),
              `inos+bodipy_pos_perc` = median(`Raw_Positive_Perc_INOS+BODIPY`),
              `inos-bodipy_pos_perc` = median(`Raw_Positive_Perc_INOS-BODIPY`),
              `bodipy_pos_perc` = median(Raw_Positive_Perc_BODIPY),
              `meanspotcount_perbodipyposcell` = median(Raw_MeanSpotCount_PerBODIPYPosCell)) %>% 
    mutate_if(is.numeric, round, 2)
  
cat(sprintf(" 
**Stimulated DMSO Raw Percentages:**\n
ARG Perc Positive: %s%% \n
INOS Perc Positive: %s%% \n
INOS+BODIPY Perc Positive: %s%% \n
INOS-BODIPY Perc Positive: %s%% \n
BODIPY Perc Positive: %s%% \n
Mean Spot Count per Bodipy Positive Cell: %.2f \n",
  dmsopercpositive$arg_pos_perc,
  dmsopercpositive$inos_pos_perc,
  dmsopercpositive$`inos+bodipy_pos_perc`,
  dmsopercpositive$`inos-bodipy_pos_perc`,
  dmsopercpositive$`bodipy_pos_perc`,
  dmsopercpositive$`meanspotcount_perbodipyposcell`))

Stimulated DMSO Raw Percentages:

ARG Perc Positive: NA%

INOS Perc Positive: NA%

INOS+BODIPY Perc Positive: NA%

INOS-BODIPY Perc Positive: NA%

BODIPY Perc Positive: NA%

Mean Spot Count per Bodipy Positive Cell: NA

ARG %-Positive Increase Hits

# set a list of variables
zlimitslist <- c( ">= 2", ">= 3", ">= 4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  
  hits_data <- ungroup(data_zscored) %>%
  filter(QCL_Sample_Number %in% hit_list & !!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>%
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Description", "Target", "Pathway", "QCL_Sample_Number", "Compound_Name", "Concentration", "Units", contains(c("Raw", "Norm", "Robust_ZScore")))

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0(feature, "$"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score >= 2

Hit dose point count

Number of unique hit compounds = 71

Hit data

Z-Score >= 3

Hit dose point count

Number of unique hit compounds = 36

Hit data

Z-Score >= 4

Hit dose point count

Number of unique hit compounds = 23

Hit data

  # cat("\n")
# }

ARG %-Positive Decrease Hits

# set a list of variables
zlimitslist <- c("<= -2", "<= -3", "<= -4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
  filter(QCL_Sample_Number %in% hit_list & !!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>%
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Compound_Name", "Concentration", "Units", "Description", "Target", "Pathway", "QCL_Sample_Number", contains(c("Raw", "Norm", "Robust_ZScore"))) 

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0(feature, "$"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score <= -2

Hit dose point count

Number of unique hit compounds = 13

Hit data

Z-Score <= -3

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -4

Hit dose point count

Number of unique hit compounds = 0

Hit data

  # cat("\n")
# }

INOS Perc Positive

Waterfall plots

Robust Z-Scores for all of the test compounds are displayed below. The plots is interactive. Hover over the points to view the compound name and other information. The number of hits (unique compound IDs) and corresponding fold change (FC) at Z-Score cut-offs of 2, 3 and 4 are printed below the plots.

# prepare the data
feat_list <- unique(data_zscored_melt$Feature)[3] #"Positive_Perc_INOS"

data_zscored_plot <- data_zscored_melt %>%
  filter(Feature == feat_list) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  unite(Treatment, c(Compound_Name, Conc), remove = FALSE) %>%
  mutate(Status = ifelse(Robust_ZScore >= 2, "Significant Increase",
                  ifelse(Robust_ZScore <= -2, "Significant Decrease", "No Change")))

# set order of status for plotting
stat_order <- c("Significant Increase", "Significant Decrease", "No Change")
data_zscored_plot$Status <- factor(data_zscored_plot$Status, levels = stat_order)

# reorder
data_zscored_plot_filt <- data_zscored_plot %>%
  group_by(Feature) %>% 
  arrange(Robust_ZScore) %>%
  mutate(Index = 1:n())

# print a plot and dataframe for each feature
# for (condition in cond_list){
  
# create the condition heading
# cat("\n")
# cat("\n### ", condition, "\n")
# cat("\n")

  for (feature in feat_list) {
  
  # create the feature heading
  cat("\n")
  cat("\n#### ", feature, "\n")
  cat("\n")

# create dot plot
snake_plot <- data_zscored_plot_filt %>%
  filter(Feature == feature) %>% 
  ggplot(., aes(x = Index, 
                y = Robust_ZScore,
                label = Compound_Name)) +
  geom_point(aes(colour = Status,
                 text = paste('Compound_Name:', Compound_Name,
                              '<br>Raw:', Raw,
                              '<br>NormToStimDMSO:', NormToStimDMSO,
                              '<br>Robust Z-Score:', Robust_ZScore,
                              '<br>PlateID:', PlateID)), 
                 alpha = 0.8, size = 1.6) +
  geom_hline(aes(yintercept = -2), linetype = "dotted", size = 0.3) +
  geom_hline(aes(yintercept = 2), linetype = "dotted", size = 0.3) +
  labs(x = "Compound",
       y = "Z-Score To Sample") +
  #scale_y_continuous(breaks = seq(-6, 6, by = 1), limits=c(-6,6)) +
  theme_bw() +
  theme(#axis.title.y = element_text(margin = margin(t = 0, r = 30, b = 0, l = 30)), 
        axis.title.y = element_text(margin = margin(t = 0, r = 10, b = 0, l = 10)),
        axis.title.x = element_text(size = 1, colour = "white"),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_text(colour = "black"),
        plot.title = element_text(margin = margin(t = 10, r = 0, b = 5, l = 0), size = 14), 
        plot.margin = margin(t = 10, r = 30, b = 5, l = 30),
        legend.title = element_blank()) +
  scale_colour_manual(values = c("#bc3754", "#2a788e", "grey61"), drop = FALSE, name = "")

print(htmltools::tagList(plotly::ggplotly(snake_plot, tooltip = "text")))

cat("\n")
cat("\n")

data_zscored_filt <- data_zscored_melt %>%
  filter(Feature == feature)

hit_count_summary_2 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_3 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 3, ">= 3",
                      ifelse(Robust_ZScore <= -3, "<= -3",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_4 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 4, ">= 4",
                      ifelse(Robust_ZScore <= -4, "<= -4",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary <- rbind(hit_count_summary_2, hit_count_summary_3, hit_count_summary_4) %>% 
  pivot_longer(c(Equivalent_FC_max, Equivalent_FC_min), 
               names_to = "Direction", values_to = "Equivalent_FC") %>% 
  filter(c(grepl('<', ZScore_Bin) & Direction == "Equivalent_FC_max") |
         c(grepl('>', ZScore_Bin) & Direction == "Equivalent_FC_min")) %>% 
  mutate(ZScore_Bin = factor(ZScore_Bin, levels = c(">= 2", ">= 3", ">= 4", "<= -2", "<= -3", "<= -4"))) %>%
  complete(ZScore_Bin, fill = list(Hit_Count = 0)) %>% 
  select(ZScore_Bin, Equivalent_FC, Hit_Count) %>%
  arrange(ZScore_Bin) %>% 
  mutate_if(is.numeric, round, 2) 
  
# display the table
print(htmltools::tagList(datatable(hit_count_summary, extensions = c('Scroller', "FixedColumns"), rownames = FALSE, 
  options = list(
  dom = 'rt',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
          pageLength = -1,
  searching = FALSE,
  scrollX = TRUE))))

  cat("\n")
  } # end the loop

Positive_Perc_INOS

# }

# Create vectors for WellID, row, and column
WellID <- character(384)
row <- character(384)
column <- integer(384)
# list_of_celllines <- unique(sort(data_zscored$Cell_Line))

# Generate WellID, row, and column values
for (i in 1:384) {
  row[i] <- LETTERS[ceiling(i/24)] # Get the row letter
  column[i] <- i %% 24              # Get the column number
  if (column[i] == 0) column[i] <- 24 # Adjust column number for multiples of 24
  WellID[i] <- paste0(row[i], sprintf("%02d", column[i])) # Combine row and column to form WellID
}

# Create a data frame
row_labels <- data.frame(WellID = WellID, row = row, column = column)

# define factor levels
counts <- c(0:9, "10+")

# function to define colour palette and matching values
scale_fill_hitcount <- function(...){
    ggplot2:::manual_scale(
        'fill',
        values = setNames(
          c("#440154", "#482878", "#3e4989", "#31688e", "#26828e",
            "#1f9e89", "#35b779", "#6ece58", "#b5de2b", "#fde725", "#f98e09"),
          counts),
        ...
    )
}

# count the plates
  plate_count <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    pull(VCFG_Plate_ID) %>%
    unique() %>%
    length()

  # count the hits per well/run
  hit_count_summary_2_well <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                        ifelse(Robust_ZScore <= -2, "<= -2",
                        "NonSig"))) %>%
    filter(!ZScore_Bin == "NonSig") %>%
    group_by(WellID, ZScore_Bin) %>%
    summarise(Hit_Count = n())

  # zscore increase
  data_plot1 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == ">= 2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = ">= 2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # zscore decrease
  data_plot2 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == "<= -2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = "<= -2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # make plot
  plot1 <- data_plot1 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    ggtitle(paste("Across", plate_count, "plate\n")) +
    labs(subtitle = " Robust_ZScore >= 2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  # make plot
  plot2 <- data_plot2 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    labs(subtitle = " Robust_ZScore <= -2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

cols <- c("#23888e", "royalblue1", "darkblue")

# count the hits per well/run
hit_count_summary_2_plate <- data_zscored_melt %>%
  filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig")),
         VCFG_Plate_ID = formatC(VCFG_Plate_ID, width = 2, format = "d", flag = "0")) %>%
  unite(VCFG_Plate_ID, c(Barcode), remove = FALSE) %>%
  group_by(VCFG_Plate_ID, Barcode, Conc, ZScore_Bin) %>%
  summarise(Hit_Count = n()) 

# create a list of plates/conc
all_plates <- distinct(hit_count_summary_2_plate[1:3])


# zscore increase
data_plot1 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == ">= 2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = ">= 2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# zscore decrease
data_plot2 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == "<= -2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = "<= -2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# make plot
plot1 <- data_plot1 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore >= 2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  # ylim(0, 15) +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_blank(),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.4, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 0, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

plot2 <- data_plot2 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore <= -2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_text(size = 11, margin = margin(t = 10, r = 0, b = 0, l = 0)),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.5, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 10, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

# } # end the loop

Hit selection INOS (Z-score)

Compounds with a negative Z score are INOS % positive increasing hits while compounds with a positive Z score are INOS % positive decreasing hits.

# include DMSO percentage positive wells
cat(sprintf(" 
**Stimulated DMSO Raw Percentages:**\n
ARG Perc Positive: %s%% \n
INOS Perc Positive: %s%% \n
INOS+BODIPY Perc Positive: %s%% \n
INOS-BODIPY Perc Positive: %s%% \n
BODIPY Perc Positive: %s%% \n
Mean Spot Count per Bodipy Positive Cell: %.2f \n",
  dmsopercpositive$arg_pos_perc,
  dmsopercpositive$inos_pos_perc,
  dmsopercpositive$`inos+bodipy_pos_perc`,
  dmsopercpositive$`inos-bodipy_pos_perc`,
  dmsopercpositive$`bodipy_pos_perc`,
  dmsopercpositive$`meanspotcount_perbodipyposcell`))

Stimulated DMSO Raw Percentages:

ARG Perc Positive: NA%

INOS Perc Positive: NA%

INOS+BODIPY Perc Positive: NA%

INOS-BODIPY Perc Positive: NA%

BODIPY Perc Positive: NA%

Mean Spot Count per Bodipy Positive Cell: NA

INOS %-Positive Increase Hits

# set a list of variables
zlimitslist <- c( ">= 2", ">= 3", ">= 4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
  filter(QCL_Sample_Number %in% hit_list & !!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>%
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Description", "Target", "Pathway", "QCL_Sample_Number", "Compound_Name", "Concentration", "Units", contains(c("Raw", "Norm", "Robust_ZScore")))

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0(feature, "$"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score >= 2

Hit dose point count

Number of unique hit compounds = 194

Hit data

Z-Score >= 3

Hit dose point count

Number of unique hit compounds = 182

Hit data

Z-Score >= 4

Hit dose point count

Number of unique hit compounds = 171

Hit data

  # cat("\n")
# }

INOS %-Positive Decrease Hits

# set a list of variables
zlimitslist <- c("<= -2", "<= -3", "<= -4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
  filter(QCL_Sample_Number %in% hit_list & !!rlang::parse_expr(paste0("Robust_ZScore_", feature, z))) %>%
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Compound_Name", "Concentration", "Units", "Description", "Target", "Pathway", "QCL_Sample_Number", contains(c("Raw", "Norm", "Robust_ZScore"))) 

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0(feature, "$"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score <= -2

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -3

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -4

Hit dose point count

Number of unique hit compounds = 0

Hit data

  # cat("\n")
# }

INOS+BODIPY Perc Positive

Waterfall plots

Robust Z-Scores for all of the test compounds are displayed below. The plots is interactive. Hover over the points to view the compound name and other information. The number of hits (unique compound IDs) and corresponding fold change (FC) at Z-Score cut-offs of 2, 3 and 4 are printed below the plots.

# prepare the data
feat_list <- unique(data_zscored_melt$Feature)[4] #"Positive_Perc_INOS+BODIPY"

data_zscored_plot <- data_zscored_melt %>%
  filter(Feature == feat_list) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  unite(Treatment, c(Compound_Name, Conc), remove = FALSE) %>%
  mutate(Status = ifelse(Robust_ZScore >= 2, "Significant Increase",
                  ifelse(Robust_ZScore <= -2, "Significant Decrease", "No Change")))

# set order of status for plotting
stat_order <- c("Significant Increase", "Significant Decrease", "No Change")
data_zscored_plot$Status <- factor(data_zscored_plot$Status, levels = stat_order)

# reorder
data_zscored_plot_filt <- data_zscored_plot %>%
  group_by(Feature) %>% 
  arrange(Robust_ZScore) %>%
  mutate(Index = 1:n())

# print a plot and dataframe for each feature
# for (condition in cond_list){
  
# create the condition heading
# cat("\n")
# cat("\n### ", condition, "\n")
# cat("\n")

  for (feature in feat_list) {
  
  # create the feature heading
  cat("\n")
  cat("\n#### ", feature, "\n")
  cat("\n")

# create dot plot
snake_plot <- data_zscored_plot_filt %>%
  filter(Feature == feature) %>% 
  ggplot(., aes(x = Index, 
                y = Robust_ZScore,
                label = Compound_Name)) +
  geom_point(aes(colour = Status,
                 text = paste('Compound_Name:', Compound_Name,
                              '<br>Raw:', Raw,
                              '<br>NormToStimDMSO:', NormToStimDMSO,
                              '<br>Robust Z-Score:', Robust_ZScore,
                              '<br>PlateID:', PlateID)), 
                 alpha = 0.8, size = 1.6) +
  geom_hline(aes(yintercept = -2), linetype = "dotted", size = 0.3) +
  geom_hline(aes(yintercept = 2), linetype = "dotted", size = 0.3) +
  labs(x = "Compound",
       y = "Z-Score To Sample") +
  #scale_y_continuous(breaks = seq(-6, 6, by = 1), limits=c(-6,6)) +
  theme_bw() +
  theme(#axis.title.y = element_text(margin = margin(t = 0, r = 30, b = 0, l = 30)), 
        axis.title.y = element_text(margin = margin(t = 0, r = 10, b = 0, l = 10)),
        axis.title.x = element_text(size = 1, colour = "white"),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_text(colour = "black"),
        plot.title = element_text(margin = margin(t = 10, r = 0, b = 5, l = 0), size = 14), 
        plot.margin = margin(t = 10, r = 30, b = 5, l = 30),
        legend.title = element_blank()) +
  scale_colour_manual(values = c("#bc3754", "#2a788e", "grey61"), drop = FALSE, name = "")

print(htmltools::tagList(plotly::ggplotly(snake_plot, tooltip = "text")))

cat("\n")
cat("\n")

data_zscored_filt <- data_zscored_melt %>%
  filter(Feature == feature)

hit_count_summary_2 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_3 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 3, ">= 3",
                      ifelse(Robust_ZScore <= -3, "<= -3",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_4 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 4, ">= 4",
                      ifelse(Robust_ZScore <= -4, "<= -4",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary <- rbind(hit_count_summary_2, hit_count_summary_3, hit_count_summary_4) %>% 
  pivot_longer(c(Equivalent_FC_max, Equivalent_FC_min), 
               names_to = "Direction", values_to = "Equivalent_FC") %>% 
  filter(c(grepl('<', ZScore_Bin) & Direction == "Equivalent_FC_max") |
         c(grepl('>', ZScore_Bin) & Direction == "Equivalent_FC_min")) %>% 
  mutate(ZScore_Bin = factor(ZScore_Bin, levels = c(">= 2", ">= 3", ">= 4", "<= -2", "<= -3", "<= -4"))) %>%
  complete(ZScore_Bin, fill = list(Hit_Count = 0)) %>% 
  select(ZScore_Bin, Equivalent_FC, Hit_Count) %>%
  arrange(ZScore_Bin) %>% 
  mutate_if(is.numeric, round, 2) 
  
# display the table
print(htmltools::tagList(datatable(hit_count_summary, extensions = c('Scroller', "FixedColumns"), rownames = FALSE, 
  options = list(
  dom = 'rt',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
          pageLength = -1,
  searching = FALSE,
  scrollX = TRUE))))

  cat("\n")
  } # end the loop

Positive_Perc_INOS+BODIPY

# }

# Create vectors for WellID, row, and column
WellID <- character(384)
row <- character(384)
column <- integer(384)
# list_of_celllines <- unique(sort(data_zscored$Cell_Line))

# Generate WellID, row, and column values
for (i in 1:384) {
  row[i] <- LETTERS[ceiling(i/24)] # Get the row letter
  column[i] <- i %% 24              # Get the column number
  if (column[i] == 0) column[i] <- 24 # Adjust column number for multiples of 24
  WellID[i] <- paste0(row[i], sprintf("%02d", column[i])) # Combine row and column to form WellID
}

# Create a data frame
row_labels <- data.frame(WellID = WellID, row = row, column = column)

# define factor levels
counts <- c(0:9, "10+")

# function to define colour palette and matching values
scale_fill_hitcount <- function(...){
    ggplot2:::manual_scale(
        'fill',
        values = setNames(
          c("#440154", "#482878", "#3e4989", "#31688e", "#26828e",
            "#1f9e89", "#35b779", "#6ece58", "#b5de2b", "#fde725", "#f98e09"),
          counts),
        ...
    )
}

# count the plates
  plate_count <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    pull(VCFG_Plate_ID) %>%
    unique() %>%
    length()

  # count the hits per well/run
  hit_count_summary_2_well <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                        ifelse(Robust_ZScore <= -2, "<= -2",
                        "NonSig"))) %>%
    filter(!ZScore_Bin == "NonSig") %>%
    group_by(WellID, ZScore_Bin) %>%
    summarise(Hit_Count = n())

  # zscore increase
  data_plot1 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == ">= 2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = ">= 2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # zscore decrease
  data_plot2 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == "<= -2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = "<= -2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # make plot
  plot1 <- data_plot1 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    ggtitle(paste("Across", plate_count, "plate\n")) +
    labs(subtitle = " Robust_ZScore >= 2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  # make plot
  plot2 <- data_plot2 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    labs(subtitle = " Robust_ZScore <= -2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

cols <- c("#23888e", "royalblue1", "darkblue")

# count the hits per well/run
hit_count_summary_2_plate <- data_zscored_melt %>%
  filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig")),
         VCFG_Plate_ID = formatC(VCFG_Plate_ID, width = 2, format = "d", flag = "0")) %>%
  unite(VCFG_Plate_ID, c(Barcode), remove = FALSE) %>%
  group_by(VCFG_Plate_ID, Barcode, Conc, ZScore_Bin) %>%
  summarise(Hit_Count = n()) 

# create a list of plates/conc
all_plates <- distinct(hit_count_summary_2_plate[1:3])


# zscore increase
data_plot1 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == ">= 2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = ">= 2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# zscore decrease
data_plot2 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == "<= -2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = "<= -2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# make plot
plot1 <- data_plot1 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore >= 2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  # ylim(0, 15) +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_blank(),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.4, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 0, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

plot2 <- data_plot2 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore <= -2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_text(size = 11, margin = margin(t = 10, r = 0, b = 0, l = 0)),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.5, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 10, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

# } # end the loop

Hit selection INOS+BODIPY (Z-score)

Compounds with a negative Z score are INOS+BODIPY % positive increasing hits while compounds with a positive Z score are INOS+BODIPY % positive decreasing hits.

# include DMSO percentage positive wells
cat(sprintf(" 
**Stimulated DMSO Raw Percentages:**\n
ARG Perc Positive: %s%% \n
INOS Perc Positive: %s%% \n
INOS+BODIPY Perc Positive: %s%% \n
INOS-BODIPY Perc Positive: %s%% \n
BODIPY Perc Positive: %s%% \n
Mean Spot Count per Bodipy Positive Cell: %.2f \n",
  dmsopercpositive$arg_pos_perc,
  dmsopercpositive$inos_pos_perc,
  dmsopercpositive$`inos+bodipy_pos_perc`,
  dmsopercpositive$`inos-bodipy_pos_perc`,
  dmsopercpositive$`bodipy_pos_perc`,
  dmsopercpositive$`meanspotcount_perbodipyposcell`))

Stimulated DMSO Raw Percentages:

ARG Perc Positive: NA%

INOS Perc Positive: NA%

INOS+BODIPY Perc Positive: NA%

INOS-BODIPY Perc Positive: NA%

BODIPY Perc Positive: NA%

Mean Spot Count per Bodipy Positive Cell: NA

INOS+BODIPY %-Positive Increase Hits

# set a list of variables
zlimitslist <- c( ">= 2", ">= 3", ">= 4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Description", "Target", "Pathway", "QCL_Sample_Number", "Compound_Name", "Concentration", "Units", contains(c("Raw", "Norm", "Robust_ZScore")))

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("Positive_Perc_INOS\\+BODIPY"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score >= 2

Hit dose point count

Number of unique hit compounds = 187

Hit data

Z-Score >= 3

Hit dose point count

Number of unique hit compounds = 176

Hit data

Z-Score >= 4

Hit dose point count

Number of unique hit compounds = 171

Hit data

  # cat("\n")
# }

INOS+BODIPY %-Positive Decrease Hits

# set a list of variables
zlimitslist <- c("<= -2", "<= -3", "<= -4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Compound_Name", "Concentration", "Units", "Description", "Target", "Pathway", "QCL_Sample_Number", contains(c("Raw", "Norm", "Robust_ZScore"))) 

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("Positive_Perc_INOS\\+BODIPY"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score <= -2

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -3

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -4

Hit dose point count

Number of unique hit compounds = 0

Hit data

  # cat("\n")
# }

INOS-BODIPY Perc Positive

Waterfall plots

Robust Z-Scores for all of the test compounds are displayed below. The plots is interactive. Hover over the points to view the compound name and other information. The number of hits (unique compound IDs) and corresponding fold change (FC) at Z-Score cut-offs of 2, 3 and 4 are printed below the plots.

# prepare the data
feat_list <- unique(data_zscored_melt$Feature)[5] #"Positive_Perc_INOS-BODIPY"

data_zscored_plot <- data_zscored_melt %>%
  filter(Feature == feat_list) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  unite(Treatment, c(Compound_Name, Conc), remove = FALSE) %>%
  mutate(Status = ifelse(Robust_ZScore >= 2, "Significant Increase",
                  ifelse(Robust_ZScore <= -2, "Significant Decrease", "No Change")))

# set order of status for plotting
stat_order <- c("Significant Increase", "Significant Decrease", "No Change")
data_zscored_plot$Status <- factor(data_zscored_plot$Status, levels = stat_order)

# reorder
data_zscored_plot_filt <- data_zscored_plot %>%
  group_by(Feature) %>% 
  arrange(Robust_ZScore) %>%
  mutate(Index = 1:n())

# print a plot and dataframe for each feature
# for (condition in cond_list){
  
# create the condition heading
# cat("\n")
# cat("\n### ", condition, "\n")
# cat("\n")

  for (feature in feat_list) {
  
  # create the feature heading
  cat("\n")
  cat("\n#### ", feature, "\n")
  cat("\n")

# create dot plot
snake_plot <- data_zscored_plot_filt %>%
  filter(Feature == feature) %>% 
  ggplot(., aes(x = Index, 
                y = Robust_ZScore,
                label = Compound_Name)) +
  geom_point(aes(colour = Status,
                 text = paste('Compound_Name:', Compound_Name,
                              '<br>Raw:', Raw,
                              '<br>NormToStimDMSO:', NormToStimDMSO,
                              '<br>Robust Z-Score:', Robust_ZScore,
                              '<br>PlateID:', PlateID)), 
                 alpha = 0.8, size = 1.6) +
  geom_hline(aes(yintercept = -2), linetype = "dotted", size = 0.3) +
  geom_hline(aes(yintercept = 2), linetype = "dotted", size = 0.3) +
  labs(x = "Compound",
       y = "Z-Score To Sample") +
  #scale_y_continuous(breaks = seq(-6, 6, by = 1), limits=c(-6,6)) +
  theme_bw() +
  theme(#axis.title.y = element_text(margin = margin(t = 0, r = 30, b = 0, l = 30)), 
        axis.title.y = element_text(margin = margin(t = 0, r = 10, b = 0, l = 10)),
        axis.title.x = element_text(size = 1, colour = "white"),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_text(colour = "black"),
        plot.title = element_text(margin = margin(t = 10, r = 0, b = 5, l = 0), size = 14), 
        plot.margin = margin(t = 10, r = 30, b = 5, l = 30),
        legend.title = element_blank()) +
  scale_colour_manual(values = c("#bc3754", "#2a788e", "grey61"), drop = FALSE, name = "")

print(htmltools::tagList(plotly::ggplotly(snake_plot, tooltip = "text")))

cat("\n")
cat("\n")

data_zscored_filt <- data_zscored_melt %>%
  filter(Feature == feature)

hit_count_summary_2 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_3 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 3, ">= 3",
                      ifelse(Robust_ZScore <= -3, "<= -3",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_4 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 4, ">= 4",
                      ifelse(Robust_ZScore <= -4, "<= -4",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary <- rbind(hit_count_summary_2, hit_count_summary_3, hit_count_summary_4) %>% 
  pivot_longer(c(Equivalent_FC_max, Equivalent_FC_min), 
               names_to = "Direction", values_to = "Equivalent_FC") %>% 
  filter(c(grepl('<', ZScore_Bin) & Direction == "Equivalent_FC_max") |
         c(grepl('>', ZScore_Bin) & Direction == "Equivalent_FC_min")) %>% 
  mutate(ZScore_Bin = factor(ZScore_Bin, levels = c(">= 2", ">= 3", ">= 4", "<= -2", "<= -3", "<= -4"))) %>%
  complete(ZScore_Bin, fill = list(Hit_Count = 0)) %>% 
  select(ZScore_Bin, Equivalent_FC, Hit_Count) %>%
  arrange(ZScore_Bin) %>% 
  mutate_if(is.numeric, round, 2) 
  
# display the table
print(htmltools::tagList(datatable(hit_count_summary, extensions = c('Scroller', "FixedColumns"), rownames = FALSE, 
  options = list(
  dom = 'rt',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
          pageLength = -1,
  searching = FALSE,
  scrollX = TRUE))))

  cat("\n")
  } # end the loop

Positive_Perc_INOS-BODIPY

# }

# Create vectors for WellID, row, and column
WellID <- character(384)
row <- character(384)
column <- integer(384)
# list_of_celllines <- unique(sort(data_zscored$Cell_Line))

# Generate WellID, row, and column values
for (i in 1:384) {
  row[i] <- LETTERS[ceiling(i/24)] # Get the row letter
  column[i] <- i %% 24              # Get the column number
  if (column[i] == 0) column[i] <- 24 # Adjust column number for multiples of 24
  WellID[i] <- paste0(row[i], sprintf("%02d", column[i])) # Combine row and column to form WellID
}

# Create a data frame
row_labels <- data.frame(WellID = WellID, row = row, column = column)

# define factor levels
counts <- c(0:9, "10+")

# function to define colour palette and matching values
scale_fill_hitcount <- function(...){
    ggplot2:::manual_scale(
        'fill',
        values = setNames(
          c("#440154", "#482878", "#3e4989", "#31688e", "#26828e",
            "#1f9e89", "#35b779", "#6ece58", "#b5de2b", "#fde725", "#f98e09"),
          counts),
        ...
    )
}

# count the plates
  plate_count <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    pull(VCFG_Plate_ID) %>%
    unique() %>%
    length()

  # count the hits per well/run
  hit_count_summary_2_well <- data_zscored_melt %>%
    filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>%
    mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                        ifelse(Robust_ZScore <= -2, "<= -2",
                        "NonSig"))) %>%
    filter(!ZScore_Bin == "NonSig") %>%
    group_by(WellID, ZScore_Bin) %>%
    summarise(Hit_Count = n())

  # zscore increase
  data_plot1 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == ">= 2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = ">= 2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # zscore decrease
  data_plot2 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == "<= -2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = "<= -2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # make plot
  plot1 <- data_plot1 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    ggtitle(paste("Across", plate_count, "plate\n")) +
    labs(subtitle = " Robust_ZScore >= 2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  # make plot
  plot2 <- data_plot2 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    labs(subtitle = " Robust_ZScore <= -2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

cols <- c("#23888e", "royalblue1", "darkblue")

# count the hits per well/run
hit_count_summary_2_plate <- data_zscored_melt %>%
  filter(Feature == feature & grepl("Unstim_M0", Barcode)) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig")),
         VCFG_Plate_ID = formatC(VCFG_Plate_ID, width = 2, format = "d", flag = "0")) %>%
  unite(VCFG_Plate_ID, c(Barcode), remove = FALSE) %>%
  group_by(VCFG_Plate_ID, Barcode, Conc, ZScore_Bin) %>%
  summarise(Hit_Count = n()) 

# create a list of plates/conc
all_plates <- distinct(hit_count_summary_2_plate[1:3])


# zscore increase
data_plot1 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == ">= 2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = ">= 2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# zscore decrease
data_plot2 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == "<= -2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = "<= -2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# make plot
plot1 <- data_plot1 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore >= 2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  # ylim(0, 15) +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_blank(),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.4, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 0, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

plot2 <- data_plot2 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore <= -2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_text(size = 11, margin = margin(t = 10, r = 0, b = 0, l = 0)),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.5, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 10, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

# } # end the loop

Hit selection INOS-BODIPY (Z-score)

Compounds with a negative Z score are INOS-BODIPY % positive increasing hits while compounds with a positive Z score are INOS-BODIPY % positive decreasing hits.

# include DMSO percentage positive wells
cat(sprintf(" 
**Stimulated DMSO Raw Percentages:**\n
ARG Perc Positive: %s%% \n
INOS Perc Positive: %s%% \n
INOS+BODIPY Perc Positive: %s%% \n
INOS-BODIPY Perc Positive: %s%% \n
BODIPY Perc Positive: %s%% \n
Mean Spot Count per Bodipy Positive Cell: %.2f \n",
  dmsopercpositive$arg_pos_perc,
  dmsopercpositive$inos_pos_perc,
  dmsopercpositive$`inos+bodipy_pos_perc`,
  dmsopercpositive$`inos-bodipy_pos_perc`,
  dmsopercpositive$`bodipy_pos_perc`,
  dmsopercpositive$`meanspotcount_perbodipyposcell`))

Stimulated DMSO Raw Percentages:

ARG Perc Positive: NA%

INOS Perc Positive: NA%

INOS+BODIPY Perc Positive: NA%

INOS-BODIPY Perc Positive: NA%

BODIPY Perc Positive: NA%

Mean Spot Count per Bodipy Positive Cell: NA

INOS-BODIPY %-Positive Increase Hits

# set a list of variables
zlimitslist <- c( ">= 2", ">= 3", ">= 4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Description", "Target", "Pathway", "QCL_Sample_Number", "Compound_Name", "Concentration", "Units", contains(c("Raw", "Norm", "Robust_ZScore")))

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("Positive_Perc_INOS-BODIPY"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score >= 2

Hit dose point count

Number of unique hit compounds = 195

Hit data

Z-Score >= 3

Hit dose point count

Number of unique hit compounds = 179

Hit data

Z-Score >= 4

Hit dose point count

Number of unique hit compounds = 167

Hit data

  # cat("\n")
# }

INOS-BODIPY %-Positive Decrease Hits

# set a list of variables
zlimitslist <- c("<= -2", "<= -3", "<= -4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Compound_Name", "Concentration", "Units", "Description", "Target", "Pathway", "QCL_Sample_Number", contains(c("Raw", "Norm", "Robust_ZScore"))) 

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("Positive_Perc_INOS-BODIPY"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score <= -2

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -3

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -4

Hit dose point count

Number of unique hit compounds = 0

Hit data

  # cat("\n")
# }

Mean Spot Count Per Cell

Waterfall plots

Robust Z-Scores for all of the test compounds are displayed below. The plots is interactive. Hover over the points to view the compound name and other information. The number of hits (unique compound IDs) and corresponding fold change (FC) at Z-Score cut-offs of 2, 3 and 4 are printed below the plots.

# prepare the data
feat_list <- unique(data_zscored_melt$Feature)[8] #"MeanSpotCount_PerCell"

data_zscored_plot <- data_zscored_melt %>%
  filter(Feature == feat_list) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  unite(Treatment, c(Compound_Name, Conc), remove = FALSE) %>%
  mutate(Status = ifelse(Robust_ZScore >= 2, "Significant Increase",
                  ifelse(Robust_ZScore <= -2, "Significant Decrease", "No Change")))

# set order of status for plotting
stat_order <- c("Significant Increase", "Significant Decrease", "No Change")
data_zscored_plot$Status <- factor(data_zscored_plot$Status, levels = stat_order)

# reorder
data_zscored_plot_filt <- data_zscored_plot %>%
  group_by(Feature) %>% 
  arrange(Robust_ZScore) %>%
  mutate(Index = 1:n())

# print a plot and dataframe for each feature
# for (condition in cond_list){
  
# create the condition heading
# cat("\n")
# cat("\n### ", condition, "\n")
# cat("\n")

  for (feature in feat_list) {
  
  # create the feature heading
  cat("\n")
  cat("\n#### ", feature, "\n")
  cat("\n")

# create dot plot
snake_plot <- data_zscored_plot_filt %>%
  filter(Feature == feature) %>% 
  ggplot(., aes(x = Index, 
                y = Robust_ZScore,
                label = Compound_Name)) +
  geom_point(aes(colour = Status,
                 text = paste('Compound_Name:', Compound_Name,
                              '<br>Raw:', Raw,
                              '<br>NormToStimDMSO:', NormToStimDMSO,
                              '<br>Robust Z-Score:', Robust_ZScore,
                              '<br>PlateID:', PlateID)), 
                 alpha = 0.8, size = 1.6) +
  geom_hline(aes(yintercept = -2), linetype = "dotted", size = 0.3) +
  geom_hline(aes(yintercept = 2), linetype = "dotted", size = 0.3) +
  labs(x = "Compound",
       y = "Z-Score To Sample") +
  #scale_y_continuous(breaks = seq(-6, 6, by = 1), limits=c(-6,6)) +
  theme_bw() +
  theme(#axis.title.y = element_text(margin = margin(t = 0, r = 30, b = 0, l = 30)), 
        axis.title.y = element_text(margin = margin(t = 0, r = 10, b = 0, l = 10)),
        axis.title.x = element_text(size = 1, colour = "white"),
        axis.text.x = element_blank(),
        axis.ticks.x = element_blank(),
        axis.text.y = element_text(colour = "black"),
        plot.title = element_text(margin = margin(t = 10, r = 0, b = 5, l = 0), size = 14), 
        plot.margin = margin(t = 10, r = 30, b = 5, l = 30),
        legend.title = element_blank()) +
  scale_colour_manual(values = c("#bc3754", "#2a788e", "grey61"), drop = FALSE, name = "")

print(htmltools::tagList(plotly::ggplotly(snake_plot, tooltip = "text")))

cat("\n")
cat("\n")

data_zscored_filt <- data_zscored_melt %>%
  filter(Feature == feature)

hit_count_summary_2 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_3 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 3, ">= 3",
                      ifelse(Robust_ZScore <= -3, "<= -3",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary_4 <- data_zscored_filt %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 4, ">= 4",
                      ifelse(Robust_ZScore <= -4, "<= -4",
                      "NonSig"))) %>% 
  filter(!ZScore_Bin == "NonSig") %>%
  group_by(ZScore_Bin) %>% 
  summarise(Equivalent_FC_max = max(NormToStimDMSO),
            Equivalent_FC_min = min(NormToStimDMSO),
            Hit_Count = n_distinct(QCL_Sample_Number))

hit_count_summary <- rbind(hit_count_summary_2, hit_count_summary_3, hit_count_summary_4) %>% 
  pivot_longer(c(Equivalent_FC_max, Equivalent_FC_min), 
               names_to = "Direction", values_to = "Equivalent_FC") %>% 
  filter(c(grepl('<', ZScore_Bin) & Direction == "Equivalent_FC_max") |
         c(grepl('>', ZScore_Bin) & Direction == "Equivalent_FC_min")) %>% 
  mutate(ZScore_Bin = factor(ZScore_Bin, levels = c(">= 2", ">= 3", ">= 4", "<= -2", "<= -3", "<= -4"))) %>%
  complete(ZScore_Bin, fill = list(Hit_Count = 0)) %>% 
  select(ZScore_Bin, Equivalent_FC, Hit_Count) %>%
  arrange(ZScore_Bin) %>% 
  mutate_if(is.numeric, round, 2) 
  
# display the table
print(htmltools::tagList(datatable(hit_count_summary, extensions = c('Scroller', "FixedColumns"), rownames = FALSE, 
  options = list(
  dom = 'rt',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
          pageLength = -1,
  searching = FALSE,
  scrollX = TRUE))))

  cat("\n")
  } # end the loop

MeanSpotCount_PerCell

# }

# Create vectors for WellID, row, and column
WellID <- character(384)
row <- character(384)
column <- integer(384)
# list_of_celllines <- unique(sort(data_zscored$Cell_Line))

# Generate WellID, row, and column values
for (i in 1:384) {
  row[i] <- LETTERS[ceiling(i/24)] # Get the row letter
  column[i] <- i %% 24              # Get the column number
  if (column[i] == 0) column[i] <- 24 # Adjust column number for multiples of 24
  WellID[i] <- paste0(row[i], sprintf("%02d", column[i])) # Combine row and column to form WellID
}

# Create a data frame
row_labels <- data.frame(WellID = WellID, row = row, column = column)

# define factor levels
counts <- c(0:9, "10+")

# function to define colour palette and matching values
scale_fill_hitcount <- function(...){
    ggplot2:::manual_scale(
        'fill',
        values = setNames(
          c("#440154", "#482878", "#3e4989", "#31688e", "#26828e",
            "#1f9e89", "#35b779", "#6ece58", "#b5de2b", "#fde725", "#f98e09"),
          counts),
        ...
    )
}

# count the plates
  plate_count <- data_zscored_melt %>%
    filter(Feature == feature) %>%
    pull(VCFG_Plate_ID) %>%
    unique() %>%
    length()

  # count the hits per well/run
  hit_count_summary_2_well <- data_zscored_melt %>%
    filter(Feature == feature) %>%
    mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                        ifelse(Robust_ZScore <= -2, "<= -2",
                        "NonSig"))) %>%
    filter(!ZScore_Bin == "NonSig") %>%
    group_by(WellID, ZScore_Bin) %>%
    summarise(Hit_Count = n())

  # zscore increase
  data_plot1 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == ">= 2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = ">= 2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # zscore decrease
  data_plot2 <- hit_count_summary_2_well %>%
    filter(ZScore_Bin == "<= -2") %>%
    left_join(row_labels, ., by = "WellID") %>%
    mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
    mutate(ZScore_Bin = "<= -2",
           Hit_Count = ifelse(Hit_Count >= 10, "10+", Hit_Count),
           Hit_Count = factor(Hit_Count, levels = counts))

  # make plot
  plot1 <- data_plot1 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    ggtitle(paste("Across", plate_count, "plate\n")) +
    labs(subtitle = " Robust_ZScore >= 2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  # make plot
  plot2 <- data_plot2 %>%
    ggplot(aes(x = factor(column), y = factor(row), fill = factor(Hit_Count, levels = counts))) +
    geom_tile(color = "black", alpha = 0.8, show.legend = TRUE) +
    scale_y_discrete(limits=rev) +
    labs(subtitle = " Robust_ZScore <= -2") +
    scale_fill_hitcount(name = "Hit Count", drop = FALSE) +
    theme(axis.title.x = element_blank(),
          axis.title.y = element_blank(),
          plot.title = element_text(size = 14),
          plot.subtitle = element_text(size = 12))

  print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

cols <- c("#23888e", "royalblue1", "darkblue")

# count the hits per well/run
hit_count_summary_2_plate <- data_zscored_melt %>%
  filter(Feature == feature) %>% 
  unite(Conc, c(Concentration, Units), sep = "") %>%
  mutate(ZScore_Bin = ifelse(Robust_ZScore >= 2, ">= 2",
                      ifelse(Robust_ZScore <= -2, "<= -2",
                      "NonSig")),
         VCFG_Plate_ID = formatC(VCFG_Plate_ID, width = 2, format = "d", flag = "0")) %>%
  unite(VCFG_Plate_ID, c(Barcode), remove = FALSE) %>%
  group_by(VCFG_Plate_ID, Barcode, Conc, ZScore_Bin) %>%
  summarise(Hit_Count = n()) 

# create a list of plates/conc
all_plates <- distinct(hit_count_summary_2_plate[1:3])


# zscore increase
data_plot1 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == ">= 2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = ">= 2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# zscore decrease
data_plot2 <- hit_count_summary_2_plate %>%
  filter(ZScore_Bin == "<= -2") %>%
  left_join(all_plates, ., by = c("Barcode", "VCFG_Plate_ID", "Conc")) %>%
  mutate_at(vars(Hit_Count), ~replace(., is.na(.), 0)) %>%
  mutate(ZScore_Bin = "<= -2",
         Conc = factor(Conc, levels = c("0.1uM", "1uM", "10uM")))

# make plot
plot1 <- data_plot1 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore >= 2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  # ylim(0, 15) +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_blank(),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.4, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 0, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

plot2 <- data_plot2 %>%
  ggplot(aes(x = factor(VCFG_Plate_ID), y = Hit_Count, fill = Conc)) +
  geom_bar(stat="identity", position=position_dodge2(reverse = TRUE), alpha = 0.7, colour = "black") +
  labs(subtitle = " Robust_ZScore <= -2",
       x = "VCFG Plate ID",
       y = "Hit Count") +
  #scale_y_continuous(breaks = seq(0, 18, by = 2), limits=c(0,NA)) +
  scale_fill_manual(values = cols, drop = FALSE) +
  guides(fill = guide_legend(reverse = TRUE)) +
  theme_bw() +
  theme(axis.title.x = element_text(size = 11, margin = margin(t = 10, r = 0, b = 0, l = 0)),
        axis.text.x = element_text(hjust = 0.8, vjust = 0.5, size = 9, angle = 45,  margin=margin(-10,0,10,0), colour = "black"),
        axis.text.y = element_text(size = 10, colour = "black"),
        axis.title.y = element_text(11, margin = margin(t = 10, r = 6, b = 0, l = 0)),
        plot.subtitle = element_text(size = 13, margin = margin(t = 10, r = 0, b = 7, l = 0)),
        legend.title = element_blank())

print((plot1 / plot2) + plot_layout(guides = 'collect', ncol = 1))

# } # end the loop

Hit selection Mean Spot Count Per Cell (Z-score)

Compounds with a negative Z score are Mean Spot Count Per Cell increasing hits while compounds with a positive Z score are Mean Spot Count Per Cell decreasing hits.

Mean Spot Count Per Cell Increase Hits

# set a list of variables
zlimitslist <- c( ">= 2", ">= 3", ">= 4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Description", "Target", "Pathway", "QCL_Sample_Number", "Compound_Name", "Concentration", "Units", contains(c("Raw", "Norm", "Robust_ZScore")))

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("MeanSpotCount_PerCell"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score >= 2

Hit dose point count

Number of unique hit compounds = 90

Hit data

Z-Score >= 3

Hit dose point count

Number of unique hit compounds = 54

Hit data

Z-Score >= 4

Hit dose point count

Number of unique hit compounds = 41

Hit data

  # cat("\n")
# }

Mean Spot Count Per Cell Decrease Hits

# set a list of variables
zlimitslist <- c("<= -2", "<= -3", "<= -4")

  for (z in zlimitslist) { # start the neg zscore loop

    # create the z score heading
  cat("\n#### Z-Score", z, "{.tabset .tabset-fade .tabset-pills}\n")

  hits_count <-
    analysed_data %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
    unite(Conc, c(Concentration, Units), sep = "", remove = FALSE) %>%
    group_by(QCL_Sample_Number, Compound_Name) %>%
    summarise(
      "0.1uM" = sum(Conc == "0.1uM"),
      "1uM" = sum(Conc == "1uM"),
      "10uM" = sum(Conc == "10uM"),
          Total = n()) %>%
    filter(Total > 0) %>%
    arrange(desc(Total))


  # create the z score heading
  cat("\n##### Hit dose point count {.tabset .tabset-fade .tabset-pills}\n")
    hit_list <- unique(hits_count$QCL_Sample_Number)
    cat("\nNumber of unique hit compounds =", length(hit_list), "\n\n")

  print(htmltools::tagList(datatable(hits_count, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
  options = list(
  dom = 'Bfltip',
  columnDefs = list(list(className = 'dt-center',targets="_all")),
  lengthMenu = list(c(5, -1), c('5', 'All')),
          pageLength = 5,
  buttons = list(list(extend = 'csv',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature)),
                 list(extend = 'excel',
                      filename = paste0(prefix, "Z", z, "_Hit_Dose_Count_", feature), title = NULL)),
  searching = FALSE,
  scrollX = TRUE))))


  cat("\n")
  cat("\n##### Hit data {.tabset .tabset-fade .tabset-pills}\n")
  cat("\n")
  hits_data <- ungroup(data_zscored) %>%
    filter(!!rlang::parse_expr(paste0("`Robust_ZScore_", feature, "`", z))) %>% # to allow Robust z score to be filter-able
  mutate(QCL_Sample_Number = factor(QCL_Sample_Number, levels = hit_list)) %>%
  arrange(QCL_Sample_Number) %>%
  # mutate(Norm_gMFI = round(Norm_gMFI, 2),
         # Robust_ZScore_gMFI = round(Robust_ZScore_gMFI, 2)) %>%
  select("Barcode", "WellID", "Compound_Name", "Concentration", "Units", "Description", "Target", "Pathway", "QCL_Sample_Number", contains(c("Raw", "Norm", "Robust_ZScore"))) 

  # display the data in a table
  print(htmltools::tagList(datatable(hits_data, extensions = c('Buttons', 'Scroller'), rownames = FALSE,
    options = list(
    dom = 'Blfrtip',
    columnDefs = list(list(className = 'dt-center',targets="_all")),
    lengthMenu = list(c(5, -1), c('5', 'All')),
            pageLength = 5,
    buttons = list(list(extend = 'csv',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature)),
                   list(extend = 'excel',
                        filename = paste0(prefix, "Z", z, "_Hit_Summary_Data_", feature), title = NULL)),
    searching = TRUE,
    scrollX = TRUE, fixedColumns = list(leftColumns = 2))) %>% 
      formatStyle(
          columns = grep(paste0("MeanSpotCount_PerCell"), names(analysed_data), value = TRUE), 
    fontWeight = "bold", backgroundColor = "lightyellow")))
  cat("\n")
  }

Z-Score <= -2

Hit dose point count

Number of unique hit compounds = 22

Hit data

Z-Score <= -3

Hit dose point count

Number of unique hit compounds = 0

Hit data

Z-Score <= -4

Hit dose point count

Number of unique hit compounds = 0

Hit data

  # cat("\n")
# }

 

---

Analysed by Ann Rann Wong

AnnRann.Wong@petermac.org

Victorian Centre for Functional Genomics

 

Session information

sessionInfo()

R version 4.4.2 (2024-10-31)
Platform: x86_64-pc-linux-gnu
Running under: Rocky Linux 9.5 (Blue Onyx)

Matrix products: default
BLAS/LAPACK: FlexiBLAS OPENBLAS-OPENMP;  LAPACK version 3.9.0

locale:
 [1] LC_CTYPE=en_AU.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_AU.UTF-8        LC_COLLATE=en_AU.UTF-8    
 [5] LC_MONETARY=en_AU.UTF-8    LC_MESSAGES=en_AU.UTF-8   
 [7] LC_PAPER=en_AU.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_AU.UTF-8 LC_IDENTIFICATION=C       

time zone: Australia/Melbourne
tzcode source: system (glibc)

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] magick_2.8.5      lubridate_1.9.4   forcats_1.0.0     stringr_1.5.1    
 [5] dplyr_1.1.4       purrr_1.0.4       readr_2.1.5       tidyr_1.3.1      
 [9] tibble_3.2.1      ggplot2_3.5.1     tidyverse_2.0.0   patchwork_1.3.0  
[13] viridis_0.6.5     viridisLite_0.4.2 reshape2_1.4.4    DT_0.33          
[17] readxl_1.4.4      data.table_1.17.0 workflowr_1.7.1  

loaded via a namespace (and not attached):
 [1] gtable_0.3.6      xfun_0.51         bslib_0.9.0       htmlwidgets_1.6.4
 [5] processx_3.8.6    callr_3.7.6       tzdb_0.4.0        vctrs_0.6.5      
 [9] tools_4.4.2       crosstalk_1.2.1   ps_1.9.0          generics_0.1.3   
[13] parallel_4.4.2    pkgconfig_2.0.3   lifecycle_1.0.4   compiler_4.4.2   
[17] farver_2.1.2      git2r_0.35.0      munsell_0.5.1     getPass_0.2-4    
[21] httpuv_1.6.15     htmltools_0.5.8.1 sass_0.4.9        lazyeval_0.2.2   
[25] yaml_2.3.10       plotly_4.10.4     crayon_1.5.3      later_1.4.1      
[29] pillar_1.10.1     jquerylib_0.1.4   whisker_0.4.1     cachem_1.1.0     
[33] mime_0.12         tidyselect_1.2.1  digest_0.6.37     stringi_1.8.4    
[37] labeling_0.4.3    rprojroot_2.0.4   fastmap_1.2.0     grid_4.4.2       
[41] colorspace_2.1-1  cli_3.6.4         magrittr_2.0.3    utf8_1.2.4       
[45] withr_3.0.2       scales_1.3.0      promises_1.3.2    bit64_4.6.0-1    
[49] timechange_0.3.0  rmarkdown_2.29    httr_1.4.7        bit_4.5.0.1      
[53] gridExtra_2.3     cellranger_1.1.0  hms_1.1.3         evaluate_1.0.3   
[57] knitr_1.49        rlang_1.1.5       Rcpp_1.0.14       glue_1.8.0       
[61] vroom_1.6.5       rstudioapi_0.17.1 jsonlite_1.9.1    R6_2.6.1         
[65] plyr_1.8.9        fs_1.6.5