diff --git a/qhtcp-workflow/apps/r/calculate_interaction_zscores.R b/qhtcp-workflow/apps/r/calculate_interaction_zscores.R index 40eede82..6962a832 100644 --- a/qhtcp-workflow/apps/r/calculate_interaction_zscores.R +++ b/qhtcp-workflow/apps/r/calculate_interaction_zscores.R @@ -150,8 +150,9 @@ load_and_filter_data <- function(easy_results_file, sd = 3) { SM = 0, OrfRep = if_else(ORF == "YDL227C", "YDL227C", OrfRep), # should these be hardcoded? conc_num = as.numeric(gsub("[^0-9\\.]", "", Conc)), - conc_num_factor = factor(as.numeric(factor(conc_num)) - 1), - conc_num_factor_num = as.numeric(conc_num_factor) + conc_num_factor_new = factor(conc_num), + conc_num_factor_zeroed = factor(as.numeric(conc_num_factor2) - 1), + conc_num_factor = as.numeric(conc_num_factor_zeroed) # for legacy purposes, neither conc_num nor a factor ) return(df) @@ -250,10 +251,10 @@ calculate_interaction_scores <- function(df, max_conc, bg_stats, Zscore_AUC = Delta_AUC / WT_sd_AUC, # Fit linear models and store in list-columns - gene_lm_L = list(lm(Delta_L ~ conc_num_factor_num, data = pick(everything()))), - gene_lm_K = list(lm(Delta_K ~ conc_num_factor_num, data = pick(everything()))), - gene_lm_r = list(lm(Delta_r ~ conc_num_factor_num, data = pick(everything()))), - gene_lm_AUC = list(lm(Delta_AUC ~ conc_num_factor_num, data = pick(everything()))), + gene_lm_L = list(lm(Delta_L ~ conc_num_factor_zeroed_num, data = pick(everything()))), + gene_lm_K = list(lm(Delta_K ~ conc_num_factor_zeroed_num, data = pick(everything()))), + gene_lm_r = list(lm(Delta_r ~ conc_num_factor_zeroed_num, data = pick(everything()))), + gene_lm_AUC = list(lm(Delta_AUC ~ conc_num_factor_zeroed_num, data = pick(everything()))), # Extract coefficients using purrr::map_dbl lm_intercept_L = map_dbl(gene_lm_L, ~ coef(.x)[1]), @@ -293,12 +294,12 @@ calculate_interaction_scores <- function(df, max_conc, bg_stats, ) calculations <- calculations %>% - mutate( - Z_lm_L = (lm_Score_L - lm_means_sds$lm_mean_L) / lm_means_sds$lm_sd_L, - Z_lm_K = (lm_Score_K - lm_means_sds$lm_mean_K) / lm_means_sds$lm_sd_K, - Z_lm_r = (lm_Score_r - lm_means_sds$lm_mean_r) / lm_means_sds$lm_sd_r, - Z_lm_AUC = (lm_Score_AUC - lm_means_sds$lm_mean_AUC) / lm_means_sds$lm_sd_AUC - ) + mutate( + Z_lm_L = (lm_Score_L - lm_means_sds$lm_mean_L) / lm_means_sds$lm_sd_L, + Z_lm_K = (lm_Score_K - lm_means_sds$lm_mean_K) / lm_means_sds$lm_sd_K, + Z_lm_r = (lm_Score_r - lm_means_sds$lm_mean_r) / lm_means_sds$lm_sd_r, + Z_lm_AUC = (lm_Score_AUC - lm_means_sds$lm_mean_AUC) / lm_means_sds$lm_sd_AUC + ) # Summarize some of the stats interactions <- calculations %>% @@ -321,7 +322,7 @@ calculate_interaction_scores <- function(df, max_conc, bg_stats, Sum_Z_Score_K = sum(Zscore_K), Sum_Z_Score_r = sum(Zscore_r), Sum_Z_Score_AUC = sum(Zscore_AUC), - + # Calculate Average Z-scores Avg_Zscore_L = Sum_Z_Score_L / num_non_removed_concs, Avg_Zscore_K = Sum_Z_Score_K / num_non_removed_concs, @@ -346,7 +347,8 @@ calculate_interaction_scores <- function(df, max_conc, bg_stats, "Exp_L", "Exp_K", "Exp_r", "Exp_AUC", "Delta_L", "Delta_K", "Delta_r", "Delta_AUC", "Zscore_L", "Zscore_K", "Zscore_r", "Zscore_AUC", - "NG", "SM", "DB") + "NG", "SM", "DB" + ) interactions <- interactions %>% select( @@ -357,30 +359,49 @@ calculate_interaction_scores <- function(df, max_conc, bg_stats, "Avg_Zscore_L", "Avg_Zscore_K", "Avg_Zscore_r", "Avg_Zscore_AUC", "lm_Score_L", "lm_Score_K", "lm_Score_r", "lm_Score_AUC", "R_Squared_L", "R_Squared_K", "R_Squared_r", "R_Squared_AUC", - "Z_lm_L", "Z_lm_K", "Z_lm_r", "Z_lm_AUC") - + "Z_lm_L", "Z_lm_K", "Z_lm_r", "Z_lm_AUC" + ) + + # Clean the original dataframe by removing overlapping columns cleaned_df <- df %>% select(-any_of( - setdiff(intersect(names(df), names(interactions)), + setdiff(intersect(names(df), names(calculations)), c("OrfRep", "Gene", "num", "conc_num", "conc_num_factor")))) - interactions_joined <- left_join(cleaned_df, interactions, by = c("OrfRep", "Gene", "num", "conc_num", "conc_num_factor")) + # Join the original dataframe with calculations + df_with_calculations <- left_join(cleaned_df, calculations, by = c("OrfRep", "Gene", "num", "conc_num", "conc_num_factor")) + + # Remove overlapping columns between df_with_calculations and interactions + df_with_calculations_clean <- df_with_calculations %>% + select(-any_of( + setdiff(intersect(names(df_with_calculations), names(interactions)), + c("OrfRep", "Gene", "num", "conc_num", "conc_num_factor")))) + + # Join with interactions to create the full dataset + full_data <- left_join(df_with_calculations_clean, interactions, by = c("OrfRep", "Gene", "num", "conc_num", "conc_num_factor")) return(list( calculations = calculations, interactions = interactions, - interactions_joined = interactions_joined)) + full_data = full_data + )) } generate_and_save_plots <- function(out_dir, filename, plot_configs) { message("Generating ", filename, ".pdf and ", filename, ".html") - # Iterate through the plot_configs (which contain both plots and grid_layout) for (config_group in plot_configs) { plot_list <- config_group$plots grid_nrow <- config_group$grid_layout$nrow grid_ncol <- config_group$grid_layout$ncol + # Set defaults if nrow or ncol are not provided + if (is.null(grid_nrow) || is.null(grid_ncol)) { + num_plots <- length(plot_list) + grid_nrow <- ifelse(is.null(grid_nrow), 1, grid_nrow) + grid_ncol <- ifelse(is.null(grid_ncol), num_plots, grid_ncol) + } + # Prepare lists to collect static and interactive plots static_plots <- list() plotly_plots <- list() @@ -419,11 +440,11 @@ generate_and_save_plots <- function(out_dir, filename, plot_configs) { # Use appropriate helper function based on plot type plot <- switch(config$plot_type, - "scatter" = generate_scatter_plot(plot, config), - "box" = generate_box_plot(plot, config), - "density" = plot + geom_density(), - "bar" = plot + geom_bar(), - plot # default case if no type matches + "scatter" = generate_scatter_plot(plot, config), + "box" = generate_box_plot(plot, config), + "density" = plot + geom_density(), + "bar" = plot + geom_bar(), + plot # default case if no type matches ) # Add title and labels @@ -462,9 +483,9 @@ generate_and_save_plots <- function(out_dir, filename, plot_configs) { # Convert to plotly object and suppress warnings here plotly_plot <- suppressWarnings({ if (length(tooltip_vars) > 0) { - ggplotly(plot, tooltip = tooltip_vars) + plotly::ggplotly(plot, tooltip = tooltip_vars) } else { - ggplotly(plot, tooltip = "none") + plotly::ggplotly(plot, tooltip = "none") } }) @@ -483,8 +504,7 @@ generate_and_save_plots <- function(out_dir, filename, plot_configs) { grid.arrange(grobs = static_plots, ncol = grid_ncol, nrow = grid_nrow) dev.off() - # Combine and save interactive HTML plot(s) - combined_plot <- subplot( + combined_plot <- plotly::subplot( plotly_plots, nrows = grid_nrow, ncols = grid_ncol, @@ -492,7 +512,8 @@ generate_and_save_plots <- function(out_dir, filename, plot_configs) { ) # Save combined HTML plot(s) - saveWidget(combined_plot, file = file.path(out_dir, paste0(filename, ".html")), selfcontained = TRUE) + html_file <- file.path(out_dir, paste0(filename, ".html")) + saveWidget(combined_plot, file = html_file, selfcontained = TRUE) } } @@ -635,8 +656,10 @@ generate_scatter_plot <- function(plot, config) { } generate_box_plot <- function(plot, config) { - plot <- plot + geom_boxplot() - + # Convert x_var to a factor within aes mapping + plot <- plot + geom_boxplot(aes(x = factor(.data[[config$x_var]]))) + + # Apply scale_x_discrete for breaks, labels, and axis label if provided if (!is.null(config$x_breaks) && !is.null(config$x_labels) && !is.null(config$x_label)) { plot <- plot + scale_x_discrete( name = config$x_label, @@ -681,7 +704,7 @@ generate_plate_analysis_plot_configs <- function(variables, stages = c("before", plot_type = plot_type, title = paste("Plate analysis by Drug Conc for", var, stage, "quality control"), error_bar = error_bar, - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", position = position, size = 0.2 ) @@ -691,75 +714,127 @@ generate_plate_analysis_plot_configs <- function(variables, stages = c("before", return(plots) } -generate_interaction_plot_configs <- function(df, limits_map = NULL, stats_df = NULL) { +generate_interaction_plot_configs <- function(df, limits_map = NULL, plot_type = "reference") { + # Define limits if not provided if (is.null(limits_map)) { limits_map <- list( L = c(0, 130), K = c(-20, 160), r = c(0, 1), - AUC = c(0, 12500) + AUC = c(0, 12500), + Delta_L = c(-60, 60), + Delta_K = c(-60, 60), + Delta_r = c(-0.6, 0.6), + Delta_AUC = c(-6000, 6000) ) } - # Ensure proper grouping by OrfRep, Gene, and num - df_filtered <- df %>% - filter( - !is.na(L) & L >= limits_map$L[1] & L <= limits_map$L[2], - !is.na(K) & K >= limits_map$K[1] & K <= limits_map$K[2], - !is.na(r) & r >= limits_map$r[1] & r <= limits_map$r[2], - !is.na(AUC) & AUC >= limits_map$AUC[1] & AUC <= limits_map$AUC[2] - ) %>% - group_by(OrfRep, Gene, num) # Group by OrfRep, Gene, and num + # Define grouping variables and filter data based on plot type + if (plot_type == "reference") { + group_vars <- c("OrfRep", "Gene", "num") + df_filtered <- df %>% + mutate( + OrfRepCombined = paste(OrfRep, Gene, num, sep = "_") + ) + } else if (plot_type == "deletion") { + group_vars <- c("OrfRep", "Gene") + df_filtered <- df %>% + mutate( + OrfRepCombined = paste(OrfRep, Gene, sep = "_") # Compare by OrfRep and Gene for deletion + ) + } - scatter_configs <- list() - box_configs <- list() + # Create a list to store all configs + configs <- list() - # Generate scatter and box plots for each variable (L, K, r, AUC) - for (var in names(limits_map)) { - scatter_configs[[length(scatter_configs) + 1]] <- list( + # Generate the first 8 scatter/box plots for L, K, r, AUC (shared between reference and deletion) + overall_vars <- c("L", "K", "r", "AUC") + for (var in overall_vars) { + y_limits <- limits_map[[var]] + + config <- list( df = df_filtered, - x_var = "conc_num", # X-axis variable - y_var = var, # Y-axis variable (Delta_L, Delta_K, Delta_r, Delta_AUC) plot_type = "scatter", + x_var = "conc_num_factor_new", + y_var = var, + x_label = unique(df_filtered$Drug)[1], title = sprintf("Scatter RF for %s with SD", var), - coord_cartesian = limits_map[[var]], # Set limits for Y-axis - annotations = list( - list(x = -0.25, y = 10, label = "NG"), - list(x = -0.25, y = 5, label = "DB"), - list(x = -0.25, y = 0, label = "SM") - ), - grid_layout = list(ncol = 4, nrow = 3) - ) - box_configs[[length(box_configs) + 1]] <- list( - df = df_filtered, - x_var = "conc_num", # X-axis variable - y_var = var, # Y-axis variable (Delta_L, Delta_K, Delta_r, Delta_AUC) - plot_type = "box", - title = sprintf("Boxplot RF for %s with SD", var), - coord_cartesian = limits_map[[var]], - grid_layout = list(ncol = 4, nrow = 3) + coord_cartesian = y_limits, + error_bar = TRUE, + x_breaks = unique(df_filtered$conc_num_factor_new), + x_labels = as.character(unique(df_filtered$conc_num)), + grid_layout = list(ncol = 2, nrow = 2) ) + configs <- append(configs, list(config)) } - # Combine scatter and box plots into grids - configs <- list( - list( - grid_layout = list(nrow = 2, ncol = 2), # Scatter plots in a 2x2 grid (for the 8 plots) - plots = scatter_configs[1:4] - ), - list( - grid_layout = list(nrow = 2, ncol = 2), # Box plots in a 2x2 grid (for the 8 plots) - plots = box_configs - ), - list( - grid_layout = list(nrow = 3, ncol = 4), # Delta_ plots in a 3x4 grid - plots = scatter_configs - ), - list( - grid_layout = list(nrow = 3, ncol = 4), # Delta_ box plots in a 3x4 grid - plots = box_configs - ) - ) + # Generate Delta comparison plots (4x3 grid for deletion and reference) + unique_groups <- df_filtered %>% select(all_of(group_vars)) %>% distinct() + + for (i in seq_len(nrow(unique_groups))) { + group <- unique_groups[i, ] + group_data <- df_filtered %>% filter(across(all_of(group_vars), ~ . == group[[cur_column()]])) + + OrfRep <- as.character(group$OrfRep) + Gene <- if ("Gene" %in% names(group)) as.character(group$Gene) else "" + num <- if ("num" %in% names(group)) as.character(group$num) else "" + OrfRepCombined <- paste(OrfRep, Gene, num, sep = "_") + + # Generate plots for Delta variables + delta_vars <- c("Delta_L", "Delta_K", "Delta_r", "Delta_AUC") + for (var in delta_vars) { + y_limits <- limits_map[[var]] + upper_y <- y_limits[2] + lower_y <- y_limits[1] + y_span <- upper_y - lower_y + + # Get WT_sd_var for error bar calculations + WT_sd_var <- paste0("WT_sd_", sub("Delta_", "", var)) + WT_sd_value <- group_data[[WT_sd_var]][1] + error_bar_ymin <- 0 - (2 * WT_sd_value) + error_bar_ymax <- 0 + (2 * WT_sd_value) + + # Set annotations (Z_Shifts, lm Z-Scores, NG, DB, SM values) + Z_Shift_var <- paste0("Z_Shift_", sub("Delta_", "", var)) + Z_lm_var <- paste0("Z_lm_", sub("Delta_", "", var)) + Z_Shift_value <- round(group_data[[Z_Shift_var]][1], 2) + Z_lm_value <- round(group_data[[Z_lm_var]][1], 2) + NG_value <- group_data$NG[1] + DB_value <- group_data$DB[1] + SM_value <- group_data$SM[1] + + annotations <- list( + list(x = 1, y = upper_y - 0.2 * y_span, label = paste("ZShift =", Z_Shift_value)), + list(x = 1, y = upper_y - 0.3 * y_span, label = paste("lm ZScore =", Z_lm_value)), + list(x = 1, y = lower_y + 0.2 * y_span, label = paste("NG =", NG_value)), + list(x = 1, y = lower_y + 0.1 * y_span, label = paste("DB =", DB_value)), + list(x = 1, y = lower_y, label = paste("SM =", SM_value)) + ) + + # Create configuration for each Delta plot + config <- list( + df = group_data, + plot_type = "scatter", + x_var = "conc_num", + y_var = var, + x_label = unique(group_data$Drug)[1], + title = paste(OrfRep, Gene, sep = " "), + coord_cartesian = y_limits, + annotations = annotations, + error_bar = TRUE, + error_bar_params = list( + ymin = error_bar_ymin, + ymax = error_bar_ymax + ), + lm_smooth = TRUE, + x_breaks = unique(group_data$conc_num_factor_new), + x_labels = as.character(unique(group_data$conc_num)), + ylim_vals = y_limits, + grid_layout = list(ncol = 4, nrow = 3) # Adjust grid layout for gene-gene comparisons + ) + configs <- append(configs, list(config)) + } + } return(configs) } @@ -826,14 +901,14 @@ generate_rank_plot_configs <- function(df, variables, is_lm = FALSE, adjust = FA annotations = list( list( x = median(df_ranked[[rank_var]], na.rm = TRUE), - y = 10, + y = max(df_ranked[[zscore_var]], na.rm = TRUE) * 0.9, label = paste("Deletion Enhancers =", num_enhancers), hjust = 0.5, vjust = 1 ), list( x = median(df_ranked[[rank_var]], na.rm = TRUE), - y = -10, + y = min(df_ranked[[zscore_var]], na.rm = TRUE) * 0.9, label = paste("Deletion Suppressors =", num_suppressors), hjust = 0.5, vjust = 0 @@ -870,7 +945,7 @@ generate_rank_plot_configs <- function(df, variables, is_lm = FALSE, adjust = FA } } - # Avg ZScore and Rank Avg ZScore Plots for r, L, K, and AUC + # Avg ZScore and Rank Avg ZScore Plots for variables for (variable in variables) { for (plot_type in c("Avg Zscore vs lm", "Rank Avg Zscore vs lm")) { @@ -894,32 +969,30 @@ generate_rank_plot_configs <- function(df, variables, is_lm = FALSE, adjust = FA # Fit the linear model lm_model <- lm(as.formula(paste(y_var, "~", x_var)), data = df_ranked) - # Extract intercept and slope from the model coefficients + # Extract intercept, slope, and R-squared from the model coefficients intercept <- coef(lm_model)[1] slope <- coef(lm_model)[2] - + r_squared <- summary(lm_model)$r.squared + + # Annotations: include R-squared in the plot + annotations <- list( + list( + x = mean(range(df_ranked[[x_var]], na.rm = TRUE)), + y = mean(range(df_ranked[[y_var]], na.rm = TRUE)), + label = paste("R-squared =", round(r_squared, 2)), + hjust = 0.5, + vjust = 1, + size = 5 + ) + ) + configs[[length(configs) + 1]] <- list( df = df_ranked, x_var = x_var, y_var = y_var, plot_type = "scatter", title = title, - annotations = list( - list( - x = median(df_ranked[[rank_var]], na.rm = TRUE), - y = 10, - label = paste("Deletion Enhancers =", num_enhancers), - hjust = 0.5, - vjust = 1 - ), - list( - x = median(df_ranked[[rank_var]], na.rm = TRUE), - y = -10, - label = paste("Deletion Suppressors =", num_suppressors), - hjust = 0.5, - vjust = 0 - ) - ), + annotations = annotations, shape = 3, size = 0.25, smooth = TRUE, @@ -936,7 +1009,7 @@ generate_rank_plot_configs <- function(df, variables, is_lm = FALSE, adjust = FA return(configs) } -generate_correlation_plot_configs <- function(df) { +generate_correlation_plot_configs <- function(df, highlight_cyan = FALSE) { # Define relationships for plotting relationships <- list( list(x = "Z_lm_L", y = "Z_lm_K", label = "Interaction L vs. Interaction K"), @@ -953,6 +1026,9 @@ generate_correlation_plot_configs <- function(df) { # Fit linear model lm_model <- lm(as.formula(paste(rel$y, "~", rel$x)), data = df) lm_summary <- summary(lm_model) + intercept <- coef(lm_model)[1] + slope <- coef(lm_model)[2] + r_squared <- lm_summary$r.squared # Construct plot configuration config <- list( @@ -965,18 +1041,18 @@ generate_correlation_plot_configs <- function(df) { y_label = paste("z-score", gsub("Z_lm_", "", rel$y)), annotations = list( list( - x = Inf, - y = Inf, - label = paste("R-squared =", round(lm_summary$r.squared, 3)), - hjust = 1.1, - vjust = 2, - size = 4, + x = mean(range(df[[rel$x]], na.rm = TRUE)), + y = mean(range(df[[rel$y]], na.rm = TRUE)), + label = paste("R-squared =", round(r_squared, 3)), + hjust = 0.5, + vjust = 1, + size = 5, color = "black" ) ), smooth = TRUE, smooth_color = "tomato3", - lm_line = list(intercept = coef(lm_model)[1], slope = coef(lm_model)[2]), + lm_line = list(intercept = intercept, slope = slope), legend_position = "right", shape = 3, size = 0.5, @@ -987,8 +1063,7 @@ generate_correlation_plot_configs <- function(df) { fill = NA, color = "grey20", alpha = 0.1 ) ), - cyan_points = TRUE, - grid_layout = list(ncol = 2, nrow = 2) + cyan_points = highlight_cyan, # Toggle cyan point highlighting ) configs[[length(configs) + 1]] <- config @@ -1023,7 +1098,7 @@ main <- function() { df_no_zeros <- df_na %>% filter(L > 0) # formerly X_noZero # Save some constants - max_conc <- max(df$conc_num_factor_num) + max_conc <- max(df$conc_num_factor_zeroed_num) l_half_median <- (median(df_above_tolerance$L, na.rm = TRUE)) / 2 k_half_median <- (median(df_above_tolerance$K, na.rm = TRUE)) / 2 @@ -1106,7 +1181,7 @@ main <- function() { plot_type = "scatter", delta_bg_point = TRUE, title = "Raw L vs K before quality control", - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", error_bar = FALSE, legend_position = "right" ) @@ -1119,7 +1194,7 @@ main <- function() { y_var = NULL, plot_type = "density", title = "Density plot for Delta Background by [Drug] (All Data)", - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", x_label = "Delta Background", y_label = "Density", error_bar = FALSE, @@ -1130,7 +1205,7 @@ main <- function() { y_var = NULL, plot_type = "bar", title = "Bar plot for Delta Background by [Drug] (All Data)", - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", x_label = "Delta Background", y_label = "Count", error_bar = FALSE, @@ -1146,7 +1221,7 @@ main <- function() { delta_bg_point = TRUE, title = paste("Raw L vs K for strains above Delta Background threshold of", round(df_above_tolerance$delta_bg_tolerance[[1]], 3), "or above"), - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", position = "jitter", annotations = list( list( @@ -1194,7 +1269,7 @@ main <- function() { plot_type = "scatter", delta_bg_point = TRUE, title = "Raw L vs K for strains falling outside 2SD of the K mean at each Conc", - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", position = "jitter", legend_position = "right" ) @@ -1208,7 +1283,7 @@ main <- function() { plot_type = "scatter", gene_point = TRUE, title = "Delta Background vs K for strains falling outside 2SD of the K mean at each Conc", - color_var = "conc_num", + color_var = "conc_num_factor_zeroed", position = "jitter", legend_position = "right" ) @@ -1305,7 +1380,7 @@ main <- function() { reference_results <- calculate_interaction_scores(df_reference_stats, max_conc, bg_stats, group_vars = c("OrfRep", "Gene", "num")) zscore_calculations_reference <- reference_results$calculations zscore_interactions_reference <- reference_results$interactions - zscore_interactions_reference_joined <- reference_results$interactions_joined + zscore_interactions_reference_joined <- reference_results$full_data message("Calculating deletion strain(s) interactions scores") df_deletion_stats <- calculate_summary_stats( @@ -1316,7 +1391,7 @@ main <- function() { deletion_results <- calculate_interaction_scores(df_deletion_stats, max_conc, bg_stats, group_vars = c("OrfRep")) zscore_calculations <- deletion_results$calculations zscore_interactions <- deletion_results$interactions - zscore_interactions_joined <- deletion_results$interactions_joined + zscore_interactions_joined <- deletion_results$full_data # Writing Z-Scores to file write.csv(zscore_calculations_reference, file = file.path(out_dir, "zscore_calculations_reference.csv"), row.names = FALSE) @@ -1326,11 +1401,11 @@ main <- function() { # Create interaction plots message("Generating reference interaction plots") - reference_plot_configs <- generate_interaction_plot_configs(zscore_interactions_reference_joined) + reference_plot_configs <- generate_interaction_plot_configs(zscore_interactions_reference_joined, plot_type = "reference") generate_and_save_plots(out_dir, "interaction_plots_reference", reference_plot_configs) message("Generating deletion interaction plots") - deletion_plot_configs <- generate_interaction_plot_configs(zscore_interactions_joined) + deletion_plot_configs <- generate_interaction_plot_configs(zscore_interactions_joined, plot_type = "deletion") generate_and_save_plots(out_dir, "interaction_plots", deletion_plot_configs) # Define conditions for enhancers and suppressors