diff --git a/.gitignore b/.gitignore index be78d5d..8d7ad8f 100644 --- a/.gitignore +++ b/.gitignore @@ -14,4 +14,7 @@ github-job-*.html job-*-logs.zip job-*-logs/ tests/testthat/Rplots.pdf +Rplots.pdf docker/project_data/*.rds +realdata-output/ +realdata-reference/ diff --git a/DESCRIPTION b/DESCRIPTION index dd9b8b2..9158491 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -55,7 +55,7 @@ Imports: Cairo (>= 1.5-12.2), circlize, clusterProfiler (>= 3.18.1), 1.6.2), stringr (>= 1.5.1), SummarizedExperiment, tibble (>= 3.1.2), tidyselect (>= 1.1.1) Suggests: BiocStyle, CALIBERrfimpute, decoupleR, dorothea, ellmer, - ExperimentHub, emmeans, ggh4x, graphlayouts, GSVA, htmlwidgets, + ExperimentHub, emmeans, ggh4x, ggrepel, graphlayouts, GSVA, htmlwidgets, knitr (>= 1.33), lme4, limma, mice, networkD3 (>= 0.4), nlme, progeny, RANN, RCy3 (>= 2.10.2), RcppAnnoy, rmarkdown, testthat (>= 3.0.0), umap, uwot, withr diff --git a/NAMESPACE b/NAMESPACE index 02c6bcd..fc914df 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -34,6 +34,7 @@ export(hc_init_save_folder) export(hc_integration) export(hc_layer_results) export(hc_list_celltype_databases) +export(hc_list_llm_models) export(hc_longitudinal_endotype_clustering) export(hc_longitudinal_meta_clustering) export(hc_longitudinal_module_cap) @@ -55,6 +56,7 @@ export(hc_pca) export(hc_pca_algo_compare) export(hc_plot_celltype_annotation_matrix) export(hc_plot_cluster_heatmap) +export(hc_plot_cluster_heatmap_view) export(hc_plot_cutoffs) export(hc_plot_deg_dist) export(hc_plot_enrichment_panels) diff --git a/NEWS.md b/NEWS.md index f98996e..b5f93f4 100644 --- a/NEWS.md +++ b/NEWS.md @@ -7,6 +7,19 @@ gap size control. - Switched functional-enrichment defaults to consistent term selection across modules and wrappers. +- Added optional DoRAG retrieval support to `hc_module_function_llm()` so LLM + module summaries can be grounded in retrieved passages and stored citations. +- Added `compare_interpretation_levels = TRUE` for side-by-side LLM + interpretations without biological context, with context, and with DoRAG. +- Added `rag_connect_timeout_sec` and `rag_continue_on_error` to make DoRAG + retrieval robust to unreachable or slow RAG servers. +- Extended `hc_plot_module_function_llm()` to plot RAG comparison fields such + as `contextual_state_rag`. +- Added common RNA-seq differential-expression packages to the Docker image, + including `DESeq2`, `limma`, `sva`, `edgeR`, and supporting visualization, + shrinkage, import, and organism annotation packages. +- Allowed `hc_split_modules()` to use one Leiden `resolution` value per + selected module. - Clarified that `hc_read_data()` now removes zero-variance genes and drops non-numeric helper columns from object-based count inputs, which can shift `hc_suggest_topvar()` inflection points slightly compared with older diff --git a/R/aaa_rng_helpers.R b/R/aaa_rng_helpers.R index b63cc58..916a017 100644 --- a/R/aaa_rng_helpers.R +++ b/R/aaa_rng_helpers.R @@ -2,7 +2,13 @@ #' @noRd .hc_with_seed <- function(seed, expr) { expr <- substitute(expr) - seed <- .hc_as_integer_safely(seed[[1]]) + # Guard against NULL / empty before indexing: `NULL[[1]]` and + # `integer(0)[[1]]` both error. A missing seed simply means "do not fix RNG". + seed <- if (is.null(seed) || length(seed) == 0L) { + NA_integer_ + } else { + .hc_as_integer_safely(seed[[1]]) + } if (length(seed) == 0L || is.na(seed)) { return(eval.parent(expr)) } diff --git a/R/hc_api.R b/R/hc_api.R index 9322b3f..f52d43f 100644 --- a/R/hc_api.R +++ b/R/hc_api.R @@ -2049,14 +2049,19 @@ hc_unsplit_modules <- function(hc, panel_key_chr <- as.character(panel_key[[1]]) title_drawn <- FALSE + is_combined_panel <- panel_key_chr %in% c("top_all_dbs", "top_all_dbs_mixed") slice_candidates <- if (identical(panel_key_chr, "top_all_dbs")) c(2L, 1L) else 1L - title_offset_mm <- max(4, 0.35 * as.numeric(fontsize)) + title_offset_mm <- if (isTRUE(is_combined_panel)) { + max(8, 0.50 * as.numeric(fontsize)) + } else { + max(4, 0.35 * as.numeric(fontsize)) + } components <- try(ComplexHeatmap::list_components(), silent = TRUE) if (inherits(components, "try-error")) { components <- character(0) } - if (identical(panel_key_chr, "top_all_dbs") && length(components) > 0) { + if (isTRUE(is_combined_panel) && length(components) > 0) { body_pattern <- paste0("^", panel_target, "_heatmap_body_[0-9]+_[0-9]+$") body_components <- components[grepl(body_pattern, components)] if (length(body_components) >= 1) { @@ -2960,6 +2965,11 @@ hc_plot_cluster_heatmap <- function(hc, file_name = "Heatmap_modules.pdf", ...) if (!inherits(hc, "HCoCenaExperiment")) { stop("`hc` must be a `HCoCenaExperiment`.") } + dot_args <- list(...) + if (!("module_label_numbering" %in% names(dot_args)) && + .hc_module_label_map_has_split_labels(hc@integration@cluster[["module_label_map"]])) { + dot_args[["module_label_numbering"]] <- "preserve_existing" + } legacy_envo <- .hc_bridge_state_env() legacy_state <- .hc_bind_bridge_hcobject( @@ -2972,7 +2982,7 @@ hc_plot_cluster_heatmap <- function(hc, file_name = "Heatmap_modules.pdf", ...) options(hcocena.suppress_legacy_warning = TRUE) on.exit(options(hcocena.suppress_legacy_warning = old_opt), add = TRUE) - base::do.call(.hc_plot_cluster_heatmap_driver, c(list(file_name = file_name), list(...))) + base::do.call(.hc_plot_cluster_heatmap_driver, c(list(file_name = file_name), dot_args)) .hc_update_hc_from_cluster_plot( hc = hc, hcobject = base::get("hcobject", envir = legacy_envo, inherits = FALSE) diff --git a/R/hc_display_helpers.R b/R/hc_display_helpers.R index aab6c3e..702c096 100644 --- a/R/hc_display_helpers.R +++ b/R/hc_display_helpers.R @@ -1,5 +1,13 @@ .hc_display_htmlwidget <- function(x) { - if (isTRUE(getOption("knitr.in.progress", FALSE)) && + display_mode <- base::tolower(base::as.character( + getOption("hcocena.htmlwidget_display", "auto") + )[[1]]) + if (!display_mode %in% c("auto", "inline", "viewer")) { + display_mode <- "auto" + } + + if (!identical(display_mode, "viewer") && + isTRUE(getOption("knitr.in.progress", FALSE)) && requireNamespace("knitr", quietly = TRUE) && isTRUE(knitr::is_html_output())) { rendered <- tryCatch( @@ -16,7 +24,42 @@ } } - utils::getFromNamespace("print.htmlwidget", "htmlwidgets")(x) + # Emit a raw text/html representation only for front-ends that consume HTML + # written to stdout (e.g. some Jupyter/IRkernel setups). This is strictly + # opt-in: RStudio notebooks and the console must NOT take this path, otherwise + # the widget is dumped as raw HTML text below the chunk instead of rendering. + emit_html <- identical(display_mode, "inline") || + isTRUE(getOption("hcocena.htmlwidget_emit_html", FALSE)) + + if (isTRUE(emit_html) && + requireNamespace("repr", quietly = TRUE)) { + displayed <- tryCatch( + { + html <- repr::repr_html(x) + if (base::is.character(html) && + base::length(html) > 0 && + base::any(base::nzchar(html))) { + base::cat(html, sep = "\n") + TRUE + } else { + FALSE + } + }, + error = function(e) FALSE + ) + if (isTRUE(displayed)) { + return(invisible(x)) + } + } + + # Default path. Dispatch through the `print` generic exactly like a bare + # `print(widget)` in a chunk would. RStudio's notebook execution renders + # widgets inline only when they go through the generic dispatch (this works + # even from inside a function); calling `print.htmlwidget` directly (e.g. via + # getFromNamespace) bypasses that hook, so the widget ends up in the Viewer + # pane instead of inline below the chunk. View defaults to interactive(), the + # same as a top-level `print(widget)`. + print(x) invisible(x) } @@ -30,12 +73,23 @@ return(base::invisible(x)) } + if (inherits(x, "hc_llm_heatmap_plot")) { + print(x) + return(base::invisible(x)) + } + if (inherits(x, "gtable")) { grid::grid.newpage() grid::grid.draw(x) return(base::invisible(x)) } + if (inherits(x, "grob")) { + grid::grid.newpage() + grid::grid.draw(x) + return(base::invisible(x)) + } + if (base::is.list(x) && !base::is.null(x$gtable) && inherits(x$gtable, "gtable")) { grid::grid.newpage() grid::grid.draw(x$gtable) diff --git a/R/hc_functional_enrichment.R b/R/hc_functional_enrichment.R index 939fa9b..5ae514e 100644 --- a/R/hc_functional_enrichment.R +++ b/R/hc_functional_enrichment.R @@ -1130,17 +1130,40 @@ functional_enrichment <- function(gene_sets = "Hallmark", ) ) } + module_label_fit <- .hc_module_label_fit_pt( + module_label_pt_size = module_label_pt_size, + module_box_width_cm = module_box_width_cm, + module_labels_display = if (heatmap_module_label_mode == "same") module_labels else NULL, + n_heat_rows = n_hc_rows, + cell_size_mm = hc_cell_mm, + module_label_fontsize = label_fontsize, + use_fontsize_request = !base::is.null(heatmap_module_label_fontsize), + fontface = "bold" + ) + module_label_pt_size_unit <- if (heatmap_module_label_mode == "same" && + base::length(module_labels) > 0 && + base::length(module_label_fit$pt_size) == base::length(module_labels)) { + grid::unit(module_label_fit$pt_size, "pt") + } else { + grid::unit(module_label_pt_size, "snpc") + } + module_label_fontsize_draw <- .hc_module_label_effective_fontsize( + module_label_fit = module_label_fit, + fallback_fontsize = label_fontsize, + fallback_pt_size = module_label_pt_size + ) module_color_map <- stats::setNames(cluster_order, cluster_order) shared_heatmap_line_lwd <- 0.5 - module_box_anno <- ComplexHeatmap::anno_simple( - cluster_order, - col = module_color_map, - pch = if (heatmap_module_label_mode == "same") module_labels else NULL, - pt_gp = grid::gpar(col = "white", fontsize = label_fontsize, fontface = "bold"), - pt_size = grid::unit(module_label_pt_size, "snpc"), - simple_anno_size = grid::unit(module_box_width_cm, "cm"), - gp = grid::gpar(col = "black", lwd = shared_heatmap_line_lwd), + module_box_anno <- .hc_module_label_box_annotation( + values = cluster_order, + colors = module_color_map, + labels = if (heatmap_module_label_mode == "same") module_labels else NULL, + label_color = "white", + label_fontsize_pt = module_label_fit$base_pt_size, + fontface = "bold", + width_cm = module_box_width_cm, + border_gp = grid::gpar(col = "black", lwd = shared_heatmap_line_lwd), which = "row" ) @@ -1806,12 +1829,19 @@ functional_enrichment <- function(gene_sets = "Hallmark", "mm", valueOnly = TRUE ) + panel_title_default_offset_mm <- 2.8 * overall_plot_scale + panel_title_top_padding_mm <- base::max( + 6 * overall_plot_scale, + panel_title_default_offset_mm + + (font_title * 25.4 / 72) + + (2.0 * overall_plot_scale) + ) draw_padding_mm <- c( - 8, - base::max(22, sig_legend_width_mm + 20), - 6, - 6 - ) * overall_plot_scale + 14 * overall_plot_scale, + 6 * overall_plot_scale, + panel_title_top_padding_mm, + base::max(22, sig_legend_width_mm + 20) * overall_plot_scale + ) draw_padding <- grid::unit(draw_padding_mm, "mm") draw_with_body_borders <- function(ht_obj, panel_title = NULL, @@ -1870,13 +1900,13 @@ functional_enrichment <- function(gene_sets = "Hallmark", invisible(drawn_ht) } pdf_height <- base::max( - 8.0, + 8.8, base::min( 24.0, ( hc_body_h_mm + base::max(20, enrichment_colname_max_cm * 10) + - 34 + 46 ) / 25.4 ) ) @@ -2249,18 +2279,32 @@ functional_enrichment <- function(gene_sets = "Hallmark", "mm", valueOnly = TRUE ) + panel_title_default_offset_mm_all <- 2.8 * overall_plot_scale + db_header_height_mm_all <- db_header_fontsize * 25.4 / 72 + panel_title_height_mm_all <- font_title * 25.4 / 72 + panel_title_header_offset_mm_all <- base::max( + panel_title_default_offset_mm_all, + db_header_height_mm_all + (4.2 * overall_plot_scale) + ) + panel_title_top_padding_mm_all <- base::max( + 6 * overall_plot_scale, + panel_title_header_offset_mm_all + + panel_title_height_mm_all + + (2.0 * overall_plot_scale) + ) draw_padding_mm_all <- c( - 8, - base::max(22, sig_legend_width_mm_all + 20), - 6, - 6 - ) * overall_plot_scale + 16 * overall_plot_scale, + 6 * overall_plot_scale, + panel_title_top_padding_mm_all, + base::max(22, sig_legend_width_mm_all + 20) * overall_plot_scale + ) draw_padding_all <- grid::unit(draw_padding_mm_all, "mm") draw_with_body_borders_all <- function(ht_obj, panel_title = NULL, panel_title_target = "enrichment_all", panel_title_slice = 1L, panel_border_slices = NULL, + panel_title_offset_mm = panel_title_default_offset_mm_all, ...) { drawn_ht <- ComplexHeatmap::draw( ht_obj, @@ -2331,12 +2375,16 @@ functional_enrichment <- function(gene_sets = "Hallmark", if (!base::is.finite(panel_title_slice) || panel_title_slice < 1) { panel_title_slice <- 1L } + panel_title_offset_mm <- .hc_first_numeric_value(panel_title_offset_mm) + if (!base::is.finite(panel_title_offset_mm) || panel_title_offset_mm <= 0) { + panel_title_offset_mm <- panel_title_default_offset_mm_all + } try( ComplexHeatmap::decorate_heatmap_body(panel_title_target, slice = panel_title_slice, { grid::grid.text( label = panel_title, x = grid::unit(0.5, "npc"), - y = grid::unit(1, "npc") + grid::unit(2.8 * overall_plot_scale, "mm"), + y = grid::unit(1, "npc") + grid::unit(panel_title_offset_mm, "mm"), just = c("center", "bottom"), gp = grid::gpar(fontsize = font_title, fontface = "bold") ) @@ -2348,13 +2396,13 @@ functional_enrichment <- function(gene_sets = "Hallmark", } pdf_height_all <- base::max( - 8.2, + 9.2, base::min( 24.0, ( hc_body_h_mm + base::max(20, enrichment_colname_max_cm_all * 10) + - 36 + 54 ) / 25.4 ) ) @@ -2398,6 +2446,7 @@ functional_enrichment <- function(gene_sets = "Hallmark", panel_title_target = "enrichment_all", panel_title_slice = all_db_title_slice, panel_border_slices = panel_border_slices_all, + panel_title_offset_mm = panel_title_header_offset_mm_all, padding = draw_padding_all ) } @@ -2408,6 +2457,7 @@ functional_enrichment <- function(gene_sets = "Hallmark", panel_title_target = "enrichment_all", panel_title_slice = all_db_title_slice, panel_border_slices = panel_border_slices_all, + panel_title_offset_mm = panel_title_header_offset_mm_all, padding = draw_padding_all ) @@ -2553,13 +2603,13 @@ functional_enrichment <- function(gene_sets = "Hallmark", enrichment_ht_all_mixed + hc_ht } pdf_height_all_mixed <- base::max( - 8.2, + 9.2, base::min( 24.0, ( hc_body_h_mm + base::max(20, enrichment_colname_max_cm_all_mixed * 10) + - 36 + 54 ) / 25.4 ) ) diff --git a/R/hc_heatmap_export_helpers.R b/R/hc_heatmap_export_helpers.R index d2220bd..a424faa 100644 --- a/R/hc_heatmap_export_helpers.R +++ b/R/hc_heatmap_export_helpers.R @@ -78,7 +78,8 @@ height, pointsize = 11, dpi = NULL, - onefile = TRUE) { + onefile = TRUE, + bg = "white") { if (requireNamespace("Cairo", quietly = TRUE)) { if (!base::is.null(dpi)) { Cairo::Cairo( @@ -88,7 +89,9 @@ pointsize = pointsize, dpi = dpi, type = "pdf", - units = "in" + units = "in", + bg = bg, + canvas = bg ) } else { Cairo::CairoPDF( @@ -96,7 +99,8 @@ width = width, height = height, pointsize = pointsize, - onefile = onefile + onefile = onefile, + bg = bg ) } } else { @@ -105,7 +109,8 @@ width = width, height = height, pointsize = pointsize, - onefile = onefile + onefile = onefile, + bg = bg ) } } @@ -127,6 +132,17 @@ ) } +.hc_draw_white_page_background <- function() { + grid::grid.rect( + x = 0.5, + y = 0.5, + width = 1, + height = 1, + gp = grid::gpar(fill = "white", col = NA) + ) + invisible(NULL) +} + .hc_export_single_page_plot <- function(file, width, height, @@ -155,6 +171,7 @@ render_page <- function(open_device) { open_device() on.exit(try(grDevices::dev.off(), silent = TRUE), add = TRUE) + .hc_draw_white_page_background() draw_fun() invisible(NULL) } @@ -190,7 +207,8 @@ height, pointsize = 11, res = 300, - draw_page_fun) { + draw_page_fun, + display = FALSE) { if (!base::is.character(file) || base::length(file) != 1 || !base::nzchar(file)) { stop("`file` must be a non-empty file path.") } @@ -218,6 +236,7 @@ ) on.exit(try(grDevices::dev.off(), silent = TRUE), add = TRUE) for (idx in base::seq_along(page_labels)) { + .hc_draw_white_page_background() draw_page_fun(idx, page_labels[[idx]]) } grDevices::dev.off() @@ -239,12 +258,24 @@ pointsize = pointsize ) on.exit(try(grDevices::dev.off(), silent = TRUE), add = TRUE) + .hc_draw_white_page_background() draw_page_fun(idx, page_labels[[idx]]) grDevices::dev.off() on.exit(NULL, add = FALSE) png_files[[idx]] <- png_file } + # Optionally replay each page on the active graphics device so the figures + # also appear inline (e.g. under an R Markdown chunk / in a notebook), in + # addition to the exported files. Gated to contexts where a display target + # exists, so batch runs do not spawn a stray Rplots.pdf. + if (isTRUE(display) && + (base::interactive() || isTRUE(base::getOption("knitr.in.progress", FALSE)))) { + for (idx in base::seq_along(page_labels)) { + draw_page_fun(idx, page_labels[[idx]]) + } + } + list( pdf = pdf_file, png = png_files @@ -271,3 +302,75 @@ } ) } + +# Save a ggplot to both a (cairo) PDF and a PNG companion in one call. +# +# Drop-in replacement for `ggplot2::ggsave(filename = "...pdf", ...)`: the PDF +# is written with `cairo_pdf` (matching the rest of the package) and a PNG of +# the same dimensions is written next to it at `res` dpi. Extra `...` arguments +# are forwarded to both saves; any `device` is ignored (PDF forces cairo_pdf, +# PNG infers from the `.png` extension). The PNG failing only warns, so a PDF is +# still produced. +.hc_ggsave_pdf_png <- function(filename, + plot, + width, + height, + units = "in", + res = 300, + ...) { + if (!base::is.character(filename) || base::length(filename) != 1 || !base::nzchar(filename)) { + stop("`filename` must be a non-empty file path.") + } + dots <- base::list(...) + dots[["device"]] <- NULL + if (!base::is.null(dots[["dpi"]])) { + res <- dots[["dpi"]] + dots[["dpi"]] <- NULL + } + + pdf_file <- .hc_export_path_with_ext(filename, "pdf") + png_file <- .hc_export_path_with_ext(filename, "png") + + pdf_args <- base::c( + base::list( + filename = pdf_file, + plot = plot, + width = width, + height = height, + units = units, + device = grDevices::cairo_pdf + ), + dots + ) + if (base::is.null(pdf_args[["bg"]])) { + pdf_args[["bg"]] <- "white" + } + base::do.call(ggplot2::ggsave, pdf_args) + + png_args <- base::c( + base::list( + filename = png_file, + plot = plot, + width = width, + height = height, + units = units, + dpi = res + ), + dots + ) + if (base::is.null(png_args[["bg"]])) { + png_args[["bg"]] <- "white" + } + tryCatch( + base::do.call(ggplot2::ggsave, png_args), + error = function(e) { + base::warning( + "Could not write PNG companion ", png_file, ": ", + base::conditionMessage(e), + call. = FALSE + ) + } + ) + + base::invisible(base::list(pdf = pdf_file, png = png_file)) +} diff --git a/R/hc_list_llm_models.R b/R/hc_list_llm_models.R new file mode 100644 index 0000000..55d6d80 --- /dev/null +++ b/R/hc_list_llm_models.R @@ -0,0 +1,218 @@ +#' List Available LLM Models +#' +#' Queries the model-listing endpoint of each selected provider and returns the +#' available models as a data frame. Supports the same providers as +#' [hc_module_function_llm()]: Gemini, Claude (Anthropic), OpenAI/ChatGPT, and a +#' local OpenAI-compatible server (vLLM). +#' +#' @param provider Character vector of providers to query. One or more of +#' `"gemini"`, `"claude"` (`"anthropic"`), `"openai"` (`"chatgpt"`), and +#' `"vllm"` (`"local"`); or `"all"` (default), which queries the three cloud +#' providers plus the local server when `base_url` is supplied. +#' @param api_key Optional API key. Only applied when a single provider is +#' requested; otherwise keys are taken from the per-provider environment +#' variables (`GEMINI_API_KEY`, `ANTHROPIC_API_KEY`, `OPENAI_API_KEY`, +#' `VLLM_API_KEY`). +#' @param base_url Base URL of the local OpenAI-compatible server. Required for +#' `"vllm"`; falls back to `VLLM_BASE_URL` or `"http://localhost:8000/v1"`. +#' @param timeout_sec Per-request timeout in seconds. Default `30`. +#' @return A data frame with columns `provider`, `model`, and `info`. Providers +#' that cannot be reached (e.g. a missing API key) are skipped with a warning +#' rather than aborting the whole call. +#' @examples +#' \dontrun{ +#' hc_list_llm_models("gemini") +#' hc_list_llm_models(c("openai", "claude")) +#' hc_list_llm_models("vllm", base_url = "http://my-server:8000/v1") +#' } +#' @export +hc_list_llm_models <- function(provider = "all", + api_key = NULL, + base_url = NULL, + timeout_sec = 30) { + providers <- base::tolower(base::as.character(provider)) + providers[providers == "chatgpt"] <- "openai" + providers[providers == "anthropic"] <- "claude" + providers[providers == "local"] <- "vllm" + if (base::any(providers == "all")) { + providers <- base::c("gemini", "claude", "openai") + if (!base::is.null(base_url) && base::nzchar(base::as.character(base_url[[1]]))) { + providers <- base::c(providers, "vllm") + } + } + providers <- base::unique(providers) + + valid <- base::c("gemini", "claude", "openai", "vllm") + unknown <- base::setdiff(providers, valid) + if (base::length(unknown) > 0) { + stop( + "Unknown provider(s): ", base::paste(unknown, collapse = ", "), + ". Use gemini, claude, openai/chatgpt, or vllm/local." + ) + } + + single <- base::length(providers) == 1L + results <- base::list() + for (p in providers) { + key_arg <- if (single) api_key else NULL + df <- tryCatch( + .hc_llm_fetch_models( + provider = p, api_key = key_arg, base_url = base_url, timeout_sec = timeout_sec + ), + error = function(e) { + warning("Skipping `", p, "`: ", base::conditionMessage(e), call. = FALSE) + NULL + } + ) + if (!base::is.null(df) && base::nrow(df) > 0) { + results[[p]] <- df + } + } + + if (base::length(results) == 0) { + return(base::data.frame( + provider = base::character(0), + model = base::character(0), + info = base::character(0), + stringsAsFactors = FALSE + )) + } + + out <- base::do.call(base::rbind, results) + out <- out[base::order(out$provider, base::tolower(out$model)), , drop = FALSE] + base::rownames(out) <- NULL + base::message( + "Found ", base::nrow(out), " model(s) across ", base::length(results), " provider(s)." + ) + out +} + + +#' First-or-default helper for nested JSON fields. +#' @noRd +.hc_llm_models_default <- function(x, default = "") { + if (base::is.null(x) || base::length(x) == 0) { + return(default) + } + base::as.character(x[[1]]) +} + + +#' Build the standard model-listing data frame. +#' @noRd +.hc_llm_models_df <- function(provider, model, info) { + base::data.frame( + provider = base::rep(provider, base::length(model)), + model = base::as.character(model), + info = base::as.character(info), + stringsAsFactors = FALSE + ) +} + + +#' GET a URL and parse the JSON body, raising informative errors on failure. +#' @noRd +.hc_llm_http_get_json <- function(url, headers = base::character(0), timeout_sec = 30) { + if (!base::requireNamespace("httr", quietly = TRUE)) { + stop("Package `httr` is required to list models.") + } + resp <- httr::GET( + url, + httr::add_headers(.headers = headers), + httr::timeout(base::max(1, base::as.numeric(timeout_sec[[1]]))) + ) + txt <- httr::content(resp, as = "text", encoding = "UTF-8") + if (httr::http_error(resp)) { + detail <- base::substr(base::gsub("\\s+", " ", base::as.character(txt)), 1, 300) + stop( + "HTTP ", httr::status_code(resp), + if (base::nzchar(detail)) base::paste0(" - ", detail) else "" + ) + } + jsonlite::fromJSON(txt, simplifyVector = FALSE) +} + + +#' Fetch the model list for a single provider. +#' @noRd +.hc_llm_fetch_models <- function(provider, api_key, base_url, timeout_sec) { + if (provider == "gemini") { + key <- .hc_llm_resolve_api_key(api_key = api_key, llm = "gemini") + js <- .hc_llm_http_get_json( + url = base::paste0( + "https://generativelanguage.googleapis.com/v1beta/models?key=", + utils::URLencode(key, reserved = TRUE) + ), + timeout_sec = timeout_sec + ) + models <- if (base::is.null(js$models)) base::list() else js$models + models <- base::Filter(function(m) { + "generateContent" %in% base::unlist(m$supportedGenerationMethods) + }, models) + if (base::length(models) == 0) { + return(NULL) + } + return(.hc_llm_models_df( + provider = "gemini", + model = base::vapply(models, function(m) base::sub("^models/", "", .hc_llm_models_default(m$name)), base::character(1)), + info = base::vapply(models, function(m) .hc_llm_models_default(m$displayName), base::character(1)) + )) + } + + if (provider == "openai") { + key <- .hc_llm_resolve_api_key(api_key = api_key, llm = "openai") + js <- .hc_llm_http_get_json( + url = "https://api.openai.com/v1/models", + headers = base::c(Authorization = base::paste("Bearer", key)), + timeout_sec = timeout_sec + ) + data <- if (base::is.null(js$data)) base::list() else js$data + if (base::length(data) == 0) { + return(NULL) + } + return(.hc_llm_models_df( + provider = "openai", + model = base::vapply(data, function(m) .hc_llm_models_default(m$id), base::character(1)), + info = base::vapply(data, function(m) .hc_llm_models_default(m$owned_by), base::character(1)) + )) + } + + if (provider == "claude") { + key <- .hc_llm_resolve_api_key(api_key = api_key, llm = "claude") + js <- .hc_llm_http_get_json( + url = "https://api.anthropic.com/v1/models?limit=100", + headers = base::c(`x-api-key` = key, `anthropic-version` = "2023-06-01"), + timeout_sec = timeout_sec + ) + data <- if (base::is.null(js$data)) base::list() else js$data + if (base::length(data) == 0) { + return(NULL) + } + return(.hc_llm_models_df( + provider = "claude", + model = base::vapply(data, function(m) .hc_llm_models_default(m$id), base::character(1)), + info = base::vapply(data, function(m) .hc_llm_models_default(m$display_name), base::character(1)) + )) + } + + if (provider == "vllm") { + burl <- .hc_llm_resolve_vllm_base_url(vllm_base_url = base_url, llm = "vllm") + key <- .hc_llm_resolve_api_key(api_key = api_key, llm = "vllm") + js <- .hc_llm_http_get_json( + url = base::paste0(burl, "/models"), + headers = base::c(Authorization = base::paste("Bearer", key)), + timeout_sec = timeout_sec + ) + data <- if (base::is.null(js$data)) base::list() else js$data + if (base::length(data) == 0) { + return(NULL) + } + return(.hc_llm_models_df( + provider = "vllm", + model = base::vapply(data, function(m) .hc_llm_models_default(m$id), base::character(1)), + info = base::vapply(data, function(m) .hc_llm_models_default(m$owned_by), base::character(1)) + )) + } + + stop("Unsupported provider: ", provider) +} diff --git a/R/hc_longitudinal_endotypes.R b/R/hc_longitudinal_endotypes.R index 6043e36..0cded9d 100644 --- a/R/hc_longitudinal_endotypes.R +++ b/R/hc_longitudinal_endotypes.R @@ -2461,7 +2461,7 @@ hc_plot_longitudinal_module_means <- function(hc, } else { as.numeric(save_height[[1]]) } - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, ".pdf")), plot = p_mod, width = w, @@ -2708,7 +2708,7 @@ hc_plot_longitudinal_module_clusters <- function(hc, } w_heat <- if (is.null(save_heatmap_width)) 13 else as.numeric(save_heatmap_width[[1]]) h_heat <- if (is.null(save_heatmap_height)) 6.5 else as.numeric(save_heatmap_height[[1]]) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Waves.pdf")), plot = p_waves, width = w_waves, @@ -2716,7 +2716,7 @@ hc_plot_longitudinal_module_clusters <- function(hc, units = "in", device = grDevices::cairo_pdf ) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Heatmap.pdf")), plot = p_heat, width = w_heat, @@ -2933,7 +2933,7 @@ hc_plot_longitudinal_endotypes <- function(hc, if (isTRUE(save_pdf)) { out_dir <- .hc_resolve_output_dir(hc) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_PCA.pdf")), plot = p_pca, width = pca_width, @@ -2941,7 +2941,7 @@ hc_plot_longitudinal_endotypes <- function(hc, units = "in", device = grDevices::cairo_pdf ) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_ModuleWaves.pdf")), plot = p_waves, width = waves_width, @@ -3636,7 +3636,7 @@ hc_plot_longitudinal_k_criterion <- function(hc, if (isTRUE(save_pdf)) { out_dir <- .hc_resolve_output_dir(hc) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Global.pdf")), plot = p_global, width = 7, @@ -3644,7 +3644,7 @@ hc_plot_longitudinal_k_criterion <- function(hc, units = "in", device = grDevices::cairo_pdf ) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Modules.pdf")), plot = p_module, width = 12, @@ -3779,7 +3779,7 @@ hc_plot_longitudinal_cap <- function(hc, if (isTRUE(save_pdf)) { out_dir <- .hc_resolve_output_dir(hc) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Heatmap.pdf")), plot = p_cap, width = 13, @@ -4072,7 +4072,7 @@ hc_plot_longitudinal_meta_embeddings <- function(hc, } if (isTRUE(save_pdf)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_PCA.pdf")), plot = p_pca, width = 7, @@ -4081,7 +4081,7 @@ hc_plot_longitudinal_meta_embeddings <- function(hc, device = grDevices::cairo_pdf ) if (!is.null(p_umap)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_UMAP.pdf")), plot = p_umap, width = 7, @@ -4091,7 +4091,7 @@ hc_plot_longitudinal_meta_embeddings <- function(hc, ) } if (!is.null(p_cross)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_CrossTab.pdf")), plot = p_cross, width = 6.5, @@ -4388,7 +4388,7 @@ hc_plot_longitudinal_meta_consensus <- function(hc, } if (isTRUE(save_pdf)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Matrix.pdf")), plot = p_heat, width = 7.5, @@ -4397,7 +4397,7 @@ hc_plot_longitudinal_meta_consensus <- function(hc, device = grDevices::cairo_pdf ) if (!is.null(p_k)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_KDiagnostics.pdf")), plot = p_k, width = 7.2, @@ -4407,7 +4407,7 @@ hc_plot_longitudinal_meta_consensus <- function(hc, ) } if (!is.null(p_stability)) { - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_Stability.pdf")), plot = p_stability, width = 6.4, @@ -4636,7 +4636,7 @@ hc_plot_longitudinal_meta_module_waves <- function(hc, } else { as.numeric(save_height[[1]]) } - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, ".pdf")), plot = p, width = w, @@ -5716,7 +5716,7 @@ hc_plot_longitudinal_enrichment_waves <- function(hc, if (isTRUE(save_pdf)) { out_dir <- .hc_resolve_output_dir(hc) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_", db_token, ".pdf")), plot = p, width = save_width, @@ -6517,7 +6517,7 @@ hc_plot_longitudinal_enrichment_meta_waves <- function(hc, if (isTRUE(save_pdf)) { out_dir <- .hc_resolve_output_dir(hc) - ggplot2::ggsave( + .hc_ggsave_pdf_png( filename = base::file.path(out_dir, base::paste0(file_prefix, "_", db_token, ".pdf")), plot = p, width = save_width, diff --git a/R/hc_longitudinal_support.R b/R/hc_longitudinal_support.R index ef66af1..2e53826 100644 --- a/R/hc_longitudinal_support.R +++ b/R/hc_longitudinal_support.R @@ -80,6 +80,20 @@ cluster_names <- base::as.character(base::unique(ids)) cluster_names <- base::setdiff(cluster_names, singletons) + + # If every cluster is a singleton there is no multi-member cluster to merge + # into. Bailing out here avoids `max(numeric(0))` (-Inf + warning) and the + # subsequent `sample(character(0), 1)` error. + if (base::length(cluster_names) == 0) { + if (base::length(singletons) > 0 && isTRUE(verbose)) { + message( + base::length(singletons), + " singletons identified but no multi-member cluster to join; left unchanged." + ) + } + return(ids) + } + connectivity <- base::numeric(base::length(cluster_names)) base::names(connectivity) <- cluster_names @@ -94,6 +108,10 @@ connectivity[[j]] <- base::mean(subSNN) } } + if (!base::any(base::is.finite(connectivity))) { + # No usable connectivity to any cluster: leave this singleton as-is. + next + } m <- base::max(connectivity, na.rm = TRUE) mi <- base::which(connectivity == m) closest_cluster <- sample(base::names(connectivity[mi]), 1) diff --git a/R/hc_module_function_gemini.R b/R/hc_module_function_gemini.R index 1ff320a..84ad126 100644 --- a/R/hc_module_function_gemini.R +++ b/R/hc_module_function_gemini.R @@ -30,6 +30,8 @@ #' @param label Optional label used for storing the result. Defaults to the #' module name or `"custom_geneset"`. Can only be used for a single module or #' a free gene set. +#' @param provider Provider-neutral alias for `llm`. When supplied it takes +#' precedence over `llm`. Accepts the same values as `llm`. #' @param llm LLM provider. Supported values are `"gemini"`, `"claude"`, #' `"anthropic"` (alias for `"claude"`), `"openai"`, `"chatgpt"` (alias for #' `"openai"`), or `"vllm"` for a local OpenAI-compatible vLLM server. @@ -45,6 +47,9 @@ #' Defaults to `"Qwen/Qwen2.5-VL-32B-Instruct"`. #' @param vllm_base_url Base URL for the local OpenAI-compatible vLLM server. #' Only used when `llm = "vllm"`. Defaults to `"http://localhost:8000/v1"`. +#' @param base_url Provider-neutral alias for the local server endpoint. When +#' supplied with `llm = "vllm"` it overrides `vllm_base_url`; for cloud +#' providers it is ignored (with a warning). #' @param model Optional generic model code. Mainly useful for OpenAI, or as a #' provider-agnostic override. #' @param max_genes Maximum number of genes sent to the API. Use `Inf` or @@ -55,6 +60,40 @@ #' requests unless overridden explicitly. #' @param pause_sec Pause in seconds between module requests. Useful for #' `module = "all"`. Use `0` or `NULL` for no pause. Default is `0`. +#' @param use_rag Logical. If `TRUE`, retrieves literature passages from the +#' DoRAG raw retrieval API and injects them into the LLM prompt as optional +#' supporting context. Default is `FALSE`. +#' @param rag_query Optional retrieval query. If `NULL`, a query is built from +#' the biological context, module label, and submitted genes. A single string +#' is reused for all modules; a named character vector can map module labels +#' to queries; a character vector with one entry per module is used in order. +#' Advanced users may pass a function with arguments `label`, `module`, +#' `genes`, and `context_text`. +#' @param rag_url DoRAG raw retrieval API endpoint. If `NULL`, uses +#' `HCOCENA_RAG_URL` or +#' `"https://limesbcnr-007901.iaas.uni-bonn.de/api/rag/query"`. +#' @param rag_category Retrieval category sent to DoRAG. Default is `"DoRAG"`. +#' @param rag_top_k Number of passages requested after reranking. Default is +#' `10`. +#' @param rag_themes Optional biological theme filter passed to DoRAG, for +#' example `"Microbiome Science"` or `c("Microbiome Science", "Innate Immunity")`. +#' @param rag_timeout_sec DoRAG request timeout in seconds. Use `0` or `NULL` +#' to disable the timeout. Default is `120`. +#' @param rag_connect_timeout_sec DoRAG connection timeout in seconds. This +#' controls how long to wait while establishing the TCP connection to the +#' RAG server. Use `0` or `NULL` to use the system default. Default is `30`. +#' @param rag_min_relevance Optional minimum rerank score. Retrieved passages +#' below this score are omitted from the prompt after retrieval. +#' @param rag_max_context_chars Maximum number of characters of formatted DoRAG +#' context injected into each LLM prompt. Use `Inf` or `NULL` to disable this +#' limit. Default is `12000`. +#' @param rag_continue_on_error Logical. If `TRUE`, a failed DoRAG request is +#' stored in the result and the LLM interpretation continues without RAG +#' passages. Default is `FALSE`. +#' @param compare_interpretation_levels Logical. If `TRUE`, stores separate +#' LLM interpretations for comparison: `without_context`, `with_context`, +#' and, when `use_rag = TRUE`, `with_rag`. The top-level `response` remains +#' the final requested interpretation. Default is `FALSE`. #' @param continue_on_error Logical. If `TRUE`, continue with the next module #' when one request fails and store the error in the result summary. #' @param save_to_hc Logical. If `TRUE`, store the result in `hc@satellite` and @@ -89,17 +128,31 @@ hc_module_function_llm <- function(hc = NULL, context = NULL, biological_context = NULL, label = NULL, + provider = NULL, llm = c("gemini", "claude", "openai", "chatgpt", "vllm"), api_key = NULL, gemini_model = NULL, claude_model = NULL, vllm_model = NULL, vllm_base_url = NULL, + base_url = NULL, model = NULL, max_genes = Inf, temperature = 0.2, timeout_sec = 60, pause_sec = 0, + use_rag = FALSE, + rag_query = NULL, + rag_url = NULL, + rag_category = "DoRAG", + rag_top_k = 10, + rag_themes = NULL, + rag_timeout_sec = 120, + rag_connect_timeout_sec = 30, + rag_min_relevance = NULL, + rag_max_context_chars = 12000, + rag_continue_on_error = FALSE, + compare_interpretation_levels = FALSE, continue_on_error = FALSE, save_to_hc = !base::is.null(hc), slot_name = "llm_module_function", @@ -117,6 +170,10 @@ hc_module_function_llm <- function(hc = NULL, stop("`hc` must be provided when `save_to_hc = TRUE`.") } + # `provider` is a sprechender alias for `llm`; when supplied it wins. + if (!base::is.null(provider) && base::nzchar(base::as.character(provider[[1]]))) { + llm <- provider + } llm <- base::tolower(base::as.character(llm[[1]])) if (!(llm %in% c("gemini", "claude", "anthropic", "openai", "chatgpt", "vllm"))) { stop("`llm` must be one of `gemini`, `claude`, `anthropic`, `openai`, `chatgpt`, or `vllm`.") @@ -180,8 +237,39 @@ hc_module_function_llm <- function(hc = NULL, if (!base::is.logical(verbose) || base::length(verbose) != 1 || base::is.na(verbose)) { stop("`verbose` must be TRUE or FALSE.") } + if (!base::is.logical(compare_interpretation_levels) || + base::length(compare_interpretation_levels) != 1 || + base::is.na(compare_interpretation_levels)) { + stop("`compare_interpretation_levels` must be TRUE or FALSE.") + } + rag_options <- .hc_llm_resolve_rag_options( + use_rag = use_rag, + rag_query = rag_query, + rag_url = rag_url, + rag_category = rag_category, + rag_top_k = rag_top_k, + rag_themes = rag_themes, + rag_timeout_sec = rag_timeout_sec, + rag_connect_timeout_sec = rag_connect_timeout_sec, + rag_min_relevance = rag_min_relevance, + rag_max_context_chars = rag_max_context_chars, + rag_continue_on_error = rag_continue_on_error + ) api_key <- .hc_llm_resolve_api_key(api_key = api_key, llm = llm) + # `base_url` is a provider-neutral alias for the local server endpoint; it + # currently maps onto the OpenAI-compatible `vllm` provider. + if (!base::is.null(base_url) && base::nzchar(base::as.character(base_url[[1]]))) { + if (llm == "vllm") { + vllm_base_url <- base_url + } else { + warning( + "`base_url` currently applies only to local `vllm` models; ignoring it for provider `", + llm, "`.", + call. = FALSE + ) + } + } vllm_base_url <- .hc_llm_resolve_vllm_base_url(vllm_base_url = vllm_base_url, llm = llm) timeout_sec <- .hc_llm_resolve_timeout( timeout_sec = timeout_sec, @@ -270,6 +358,8 @@ hc_module_function_llm <- function(hc = NULL, continue_on_error = continue_on_error, system_instruction = system_instruction, response_schema = .hc_llm_response_schema(), + rag_options = rag_options, + compare_interpretation_levels = compare_interpretation_levels, verbose = verbose ) @@ -330,6 +420,8 @@ hc_module_function_vllm <- function(...) { continue_on_error, system_instruction, response_schema, + rag_options, + compare_interpretation_levels, verbose) { results <- vector("list", length = base::length(gene_infos)) for (i in base::seq_along(gene_infos)) { @@ -353,6 +445,10 @@ hc_module_function_vllm <- function(...) { timeout_sec = timeout_sec, system_instruction = system_instruction, response_schema = response_schema, + rag_options = rag_options, + compare_interpretation_levels = compare_interpretation_levels, + index = i, + n_inputs = base::length(gene_infos), verbose = verbose ), error = function(e) { @@ -531,6 +627,10 @@ hc_module_function_vllm <- function(...) { timeout_sec, system_instruction, response_schema, + rag_options = NULL, + compare_interpretation_levels = FALSE, + index = 1L, + n_inputs = 1L, verbose) { genes_all <- .hc_gemini_normalize_genes(gene_info$genes) if (base::length(genes_all) == 0) { @@ -539,6 +639,10 @@ hc_module_function_vllm <- function(...) { genes_use <- .hc_llm_limit_genes(genes = genes_all, max_genes = max_genes) truncated <- base::length(genes_use) < base::length(genes_all) + rag_result <- NULL + rag_query_used <- NULL + rag_context_text <- NULL + rag_error_message <- NA_character_ if (isTRUE(verbose)) { if (truncated) { @@ -559,16 +663,238 @@ hc_module_function_vllm <- function(...) { } } + if (!base::is.null(rag_options)) { + rag_query_used <- .hc_llm_resolve_rag_query( + rag_query = rag_options$query, + label = label, + module = gene_info$module, + genes = genes_use, + context_text = context_text, + index = index, + n_inputs = n_inputs + ) + if (isTRUE(verbose)) { + message( + "DoRAG retrieval: requesting ", + rag_options$top_k, + " passages for `", label, "`." + ) + } + rag_raw <- tryCatch( + .hc_llm_request_rag( + query = rag_query_used, + url = rag_options$url, + category = rag_options$category, + top_k = rag_options$top_k, + themes = rag_options$themes, + timeout_sec = rag_options$timeout_sec, + connect_timeout_sec = rag_options$connect_timeout_sec + ), + error = function(e) { + if (!isTRUE(rag_options$continue_on_error)) { + stop(e) + } + rag_error_message <<- base::conditionMessage(e) + if (isTRUE(verbose)) { + message( + "DoRAG retrieval failed for `", + label, + "`; continuing without RAG passages: ", + rag_error_message + ) + } + NULL + } + ) + if (base::is.null(rag_raw)) { + rag_result <- list( + query = rag_query_used, + answer = "", + context = list(), + cited_papers = list(), + status = "error", + error_message = rag_error_message + ) + } else { + rag_result <- .hc_llm_prepare_rag_result( + rag_result = rag_raw, + min_relevance = rag_options$min_relevance + ) + rag_context_text <- .hc_llm_format_rag_context( + rag_result = rag_result, + max_chars = rag_options$max_context_chars + ) + rag_result$status <- "ok" + rag_result$error_message <- NA_character_ + } + if (isTRUE(verbose)) { + message( + "DoRAG retrieval: using ", + .hc_llm_rag_context_count(rag_result), + " retrieved passages for `", label, "`." + ) + } + } + + run_level <- function(level, level_context_text, level_rag_context_text) { + .hc_llm_interpretation_level( + level = level, + label = label, + genes = genes_use, + total_gene_count = base::length(genes_all), + context_text = level_context_text, + truncated = truncated, + rag_context_text = level_rag_context_text, + llm = llm, + api_key = api_key, + model = model, + vllm_base_url = vllm_base_url, + temperature = temperature, + timeout_sec = timeout_sec, + system_instruction = system_instruction, + response_schema = response_schema + ) + } + + has_rag_context <- !base::is.null(rag_context_text) && + base::nzchar(base::as.character(rag_context_text[[1]])) + primary_level <- if (has_rag_context) { + "with_rag" + } else if (base::nzchar(context_text)) { + "with_context" + } else { + "without_context" + } + + interpretation_levels <- NULL + if (isTRUE(compare_interpretation_levels)) { + if (isTRUE(verbose)) { + message("LLM module summary: running comparison levels for `", label, "`.") + } + interpretation_levels <- list( + without_context = run_level( + level = "without_context", + level_context_text = "", + level_rag_context_text = NULL + ), + with_context = run_level( + level = "with_context", + level_context_text = context_text, + level_rag_context_text = NULL + ) + ) + if (has_rag_context) { + interpretation_levels$with_rag <- run_level( + level = "with_rag", + level_context_text = context_text, + level_rag_context_text = rag_context_text + ) + } + primary <- interpretation_levels[[primary_level]] + } else { + primary <- run_level( + level = primary_level, + level_context_text = context_text, + level_rag_context_text = rag_context_text + ) + } + + out <- list( + label = label, + module = gene_info$module, + genes_input = genes_all, + genes_sent = genes_use, + gene_count_input = base::length(genes_all), + gene_count_sent = base::length(genes_use), + truncated = truncated, + context = if (base::nzchar(context_text)) context_text else NULL, + llm = llm, + model = model, + status = "ok", + error_message = NA_character_, + prompt = primary$prompt, + response = primary$response, + response_text = primary$response_text, + raw_response_text = primary$raw_response_text, + rag = rag_result, + rag_query = rag_query_used, + rag_error_message = rag_error_message, + rag_context_text = if (!base::is.null(rag_context_text) && base::nzchar(rag_context_text)) rag_context_text else NULL, + primary_interpretation_level = primary_level, + timestamp = base::as.character(Sys.time()) + ) + if (!base::is.null(interpretation_levels)) { + out$interpretation_levels <- interpretation_levels + } + out +} + +.hc_llm_interpretation_level <- function(level, + label, + genes, + total_gene_count, + context_text, + truncated, + rag_context_text, + llm, + api_key, + model, + vllm_base_url, + temperature, + timeout_sec, + system_instruction, + response_schema) { prompt <- .hc_gemini_build_prompt( label = label, - genes = genes_use, - total_gene_count = base::length(genes_all), + genes = genes, + total_gene_count = total_gene_count, context_text = context_text, truncated = truncated, + rag_context_text = rag_context_text, llm = llm ) - req_result <- if (llm == "gemini") { + req_result <- .hc_llm_request_by_provider( + llm = llm, + api_key = api_key, + model = model, + prompt = prompt, + system_instruction = system_instruction, + response_schema = response_schema, + temperature = temperature, + timeout_sec = timeout_sec, + vllm_base_url = vllm_base_url + ) + result_text <- req_result$result_text + + list( + level = level, + context = if (base::nzchar(context_text)) context_text else NULL, + rag_used = !base::is.null(rag_context_text) && + base::nzchar(base::as.character(rag_context_text[[1]])), + rag_context_text = if (!base::is.null(rag_context_text) && + base::nzchar(base::as.character(rag_context_text[[1]]))) { + base::as.character(rag_context_text[[1]]) + } else { + NULL + }, + prompt = prompt, + response = .hc_llm_parse_module_response(result_text), + response_text = result_text, + raw_response_text = req_result$raw_response_text + ) +} + +.hc_llm_request_by_provider <- function(llm, + api_key, + model, + prompt, + system_instruction, + response_schema, + temperature, + timeout_sec, + vllm_base_url) { + if (llm == "gemini") { .hc_llm_request_gemini( api_key = api_key, model = model, @@ -608,9 +934,10 @@ hc_module_function_vllm <- function(...) { timeout_sec = timeout_sec ) } +} - result_text <- req_result$result_text - parsed_result <- tryCatch( +.hc_llm_parse_module_response <- function(result_text) { + tryCatch( jsonlite::fromJSON(result_text, simplifyVector = TRUE), error = function(e) { list( @@ -621,26 +948,6 @@ hc_module_function_vllm <- function(...) { ) } ) - - list( - label = label, - module = gene_info$module, - genes_input = genes_all, - genes_sent = genes_use, - gene_count_input = base::length(genes_all), - gene_count_sent = base::length(genes_use), - truncated = truncated, - context = if (base::nzchar(context_text)) context_text else NULL, - llm = llm, - model = model, - status = "ok", - error_message = NA_character_, - prompt = prompt, - response = parsed_result, - response_text = result_text, - raw_response_text = req_result$raw_response_text, - timestamp = base::as.character(Sys.time()) - ) } .hc_llm_limit_genes <- function(genes, max_genes) { @@ -650,6 +957,512 @@ hc_module_function_vllm <- function(...) { utils::head(genes, base::as.integer(max_genes[[1]])) } +.hc_llm_default_rag_url <- function() { + Sys.getenv("HCOCENA_RAG_URL", unset = "https://limesbcnr-007901.iaas.uni-bonn.de/api/rag/query") +} + +.hc_llm_resolve_rag_options <- function(use_rag, + rag_query, + rag_url, + rag_category, + rag_top_k, + rag_themes, + rag_timeout_sec, + rag_connect_timeout_sec, + rag_min_relevance, + rag_max_context_chars, + rag_continue_on_error) { + if (!base::is.logical(use_rag) || base::length(use_rag) != 1 || base::is.na(use_rag)) { + stop("`use_rag` must be TRUE or FALSE.") + } + if (!isTRUE(use_rag)) { + return(NULL) + } + + if (base::is.null(rag_url) || !base::nzchar(base::as.character(rag_url[[1]]))) { + rag_url <- .hc_llm_default_rag_url() + } else { + rag_url <- base::as.character(rag_url[[1]]) + } + rag_url <- base::trimws(rag_url) + if (!base::nzchar(rag_url)) { + stop("`rag_url` must be a non-empty URL when `use_rag = TRUE`.") + } + + if (base::is.null(rag_category) || !base::nzchar(base::as.character(rag_category[[1]]))) { + stop("`rag_category` must be a non-empty character scalar.") + } + rag_category <- base::trimws(base::as.character(rag_category[[1]])) + if (!base::nzchar(rag_category)) { + stop("`rag_category` must be a non-empty character scalar.") + } + + if (!base::is.numeric(rag_top_k) || base::length(rag_top_k) != 1 || + base::is.na(rag_top_k[[1]]) || !base::is.finite(rag_top_k[[1]]) || + rag_top_k[[1]] <= 0) { + stop("`rag_top_k` must be a single positive number.") + } + rag_top_k <- base::as.integer(rag_top_k[[1]]) + + if (!base::is.null(rag_timeout_sec)) { + if (!base::is.numeric(rag_timeout_sec) || base::length(rag_timeout_sec) != 1 || + base::is.na(rag_timeout_sec[[1]]) || !base::is.finite(rag_timeout_sec[[1]])) { + stop("`rag_timeout_sec` must be NULL or a single finite numeric value >= 0.") + } + } + rag_timeout_sec <- .hc_llm_normalize_timeout(rag_timeout_sec) + + if (!base::is.null(rag_connect_timeout_sec)) { + if (!base::is.numeric(rag_connect_timeout_sec) || base::length(rag_connect_timeout_sec) != 1 || + base::is.na(rag_connect_timeout_sec[[1]]) || !base::is.finite(rag_connect_timeout_sec[[1]])) { + stop("`rag_connect_timeout_sec` must be NULL or a single finite numeric value >= 0.") + } + } + rag_connect_timeout_sec <- .hc_llm_normalize_timeout(rag_connect_timeout_sec) + + if (!base::is.null(rag_min_relevance)) { + if (!base::is.numeric(rag_min_relevance) || base::length(rag_min_relevance) != 1 || + base::is.na(rag_min_relevance[[1]]) || !base::is.finite(rag_min_relevance[[1]])) { + stop("`rag_min_relevance` must be NULL or a single finite numeric value.") + } + rag_min_relevance <- base::as.numeric(rag_min_relevance[[1]]) + } + + if (!base::is.null(rag_max_context_chars)) { + if (!base::is.numeric(rag_max_context_chars) || base::length(rag_max_context_chars) != 1 || + base::is.na(rag_max_context_chars[[1]]) || rag_max_context_chars[[1]] <= 0) { + stop("`rag_max_context_chars` must be NULL, Inf, or a single positive number.") + } + if (!base::is.finite(rag_max_context_chars[[1]])) { + rag_max_context_chars <- NULL + } else { + rag_max_context_chars <- base::as.integer(rag_max_context_chars[[1]]) + } + } + + if (!base::is.null(rag_themes)) { + rag_themes <- base::as.character(rag_themes) + rag_themes <- base::trimws(rag_themes) + rag_themes <- base::unique(rag_themes[!base::is.na(rag_themes) & base::nzchar(rag_themes)]) + if (base::length(rag_themes) == 0) { + rag_themes <- NULL + } + } + + if (!base::is.logical(rag_continue_on_error) || + base::length(rag_continue_on_error) != 1 || + base::is.na(rag_continue_on_error)) { + stop("`rag_continue_on_error` must be TRUE or FALSE.") + } + + if (!base::is.null(rag_query) && !base::is.function(rag_query)) { + query_names <- base::names(rag_query) + rag_query <- base::as.character(rag_query) + rag_query <- base::trimws(rag_query) + keep <- !base::is.na(rag_query) & base::nzchar(rag_query) + if (base::length(query_names) == base::length(rag_query)) { + query_names <- base::as.character(query_names) + query_names[base::is.na(query_names)] <- "" + query_names <- query_names[keep] + } else { + query_names <- NULL + } + rag_query <- rag_query[keep] + if (base::length(rag_query) == 0) { + stop("`rag_query` must contain at least one non-empty query when provided.") + } + if (!base::is.null(query_names)) { + base::names(rag_query) <- query_names + } + } + + list( + query = rag_query, + url = rag_url, + category = rag_category, + top_k = rag_top_k, + themes = rag_themes, + timeout_sec = rag_timeout_sec, + connect_timeout_sec = rag_connect_timeout_sec, + min_relevance = rag_min_relevance, + max_context_chars = rag_max_context_chars, + continue_on_error = rag_continue_on_error + ) +} + +.hc_llm_resolve_rag_query <- function(rag_query, + label, + module, + genes, + context_text, + index, + n_inputs) { + if (base::is.null(rag_query)) { + return(.hc_llm_auto_rag_query(label = label, genes = genes, context_text = context_text)) + } + + if (base::is.function(rag_query)) { + query <- rag_query( + label = label, + module = module, + genes = genes, + context_text = context_text + ) + return(.hc_llm_normalize_rag_query_scalar(query, context = "`rag_query` function result")) + } + + query_names <- base::names(rag_query) + if (!base::is.null(query_names) && base::length(query_names) == base::length(rag_query)) { + query_names <- base::as.character(query_names) + idx <- base::match(base::as.character(label[[1]]), query_names) + if (base::is.na(idx) && !base::is.null(module)) { + idx <- base::match(base::as.character(module[[1]]), query_names) + } + if (!base::is.na(idx)) { + return(.hc_llm_normalize_rag_query_scalar(rag_query[[idx]], context = "`rag_query`")) + } + } + + if (base::length(rag_query) == 1) { + return(.hc_llm_normalize_rag_query_scalar(rag_query[[1]], context = "`rag_query`")) + } + + if (base::length(rag_query) == n_inputs) { + return(.hc_llm_normalize_rag_query_scalar(rag_query[[index]], context = "`rag_query`")) + } + + stop( + "`rag_query` must be NULL, a single string, a named vector matching module labels, ", + "or a vector with one query per module.", + call. = FALSE + ) +} + +.hc_llm_auto_rag_query <- function(label, genes, context_text, max_genes = 60) { + genes <- .hc_gemini_normalize_genes(genes) + gene_text <- if (base::length(genes) > 0) { + base::paste(utils::head(genes, max_genes), collapse = ", ") + } else { + "" + } + parts <- base::c( + "Transcriptomic module biological function", + if (base::nzchar(context_text)) base::paste0("Biological context: ", context_text) else NULL, + if (!base::is.null(label) && base::nzchar(base::as.character(label[[1]]))) base::paste0("Module label: ", base::as.character(label[[1]])) else NULL, + if (base::nzchar(gene_text)) base::paste0("Genes: ", gene_text) else NULL + ) + stringr::str_squish(base::paste(parts, collapse = "\n")) +} + +.hc_llm_normalize_rag_query_scalar <- function(query, context = "`rag_query`") { + if (base::is.null(query) || base::length(query) == 0) { + stop(context, " must resolve to a non-empty character scalar.", call. = FALSE) + } + query <- base::as.character(query[[1]]) + query <- base::trimws(query) + if (base::is.na(query) || !base::nzchar(query)) { + stop(context, " must resolve to a non-empty character scalar.", call. = FALSE) + } + query +} + +.hc_llm_request_rag <- function(query, + url, + category, + top_k, + themes, + timeout_sec, + connect_timeout_sec = NULL) { + if (!base::requireNamespace("httr", quietly = TRUE)) { + stop("Package `httr` is required for DoRAG retrieval.", call. = FALSE) + } + + payload <- list( + query = query, + category = category, + top_k = base::as.integer(top_k) + ) + if (!base::is.null(themes) && base::length(themes) > 0) { + payload$themes <- themes + } + + configs <- list(httr::add_headers("Content-Type" = "application/json")) + if (!base::is.null(timeout_sec)) { + configs <- base::c(configs, list(httr::timeout(base::as.numeric(timeout_sec[[1]])))) + } + if (!base::is.null(connect_timeout_sec)) { + configs <- base::c(configs, list(httr::config(connecttimeout = base::as.numeric(connect_timeout_sec[[1]])))) + } + + resp <- tryCatch( + do.call( + httr::POST, + base::c( + list( + url = url, + body = payload, + encode = "json" + ), + configs + ) + ), + error = function(e) { + stop("DoRAG retrieval request failed: ", base::conditionMessage(e), call. = FALSE) + } + ) + + txt <- tryCatch( + httr::content(resp, as = "text", encoding = "UTF-8"), + error = function(e) { + stop("Could not read DoRAG retrieval response: ", base::conditionMessage(e), call. = FALSE) + } + ) + txt <- base::paste(base::as.character(txt), collapse = "") + + if (httr::http_error(resp)) { + stop( + "DoRAG retrieval failed with HTTP ", + httr::status_code(resp), + if (base::nzchar(txt)) base::paste0(": ", stringr::str_trunc(stringr::str_squish(txt), 300)) else "", + call. = FALSE + ) + } + if (!base::nzchar(txt)) { + stop("DoRAG retrieval response was empty.", call. = FALSE) + } + + out <- tryCatch( + jsonlite::fromJSON(txt, simplifyVector = FALSE), + error = function(e) { + stop("Could not parse DoRAG retrieval JSON: ", base::conditionMessage(e), call. = FALSE) + } + ) + if (!base::is.list(out)) { + stop("DoRAG retrieval response did not contain a JSON object.", call. = FALSE) + } + out$raw_response_text <- txt + out +} + +.hc_llm_prepare_rag_result <- function(rag_result, min_relevance = NULL) { + if (base::is.null(rag_result) || !base::is.list(rag_result)) { + return(NULL) + } + + context_entries <- .hc_llm_rag_context_list(rag_result[["context"]]) + context_entries <- lapply(context_entries, .hc_llm_normalize_rag_context_entry) + context_entries <- context_entries[base::vapply( + context_entries, + function(x) base::nzchar(x$chunk), + FUN.VALUE = base::logical(1) + )] + + if (!base::is.null(min_relevance)) { + context_entries <- context_entries[base::vapply( + context_entries, + function(x) !base::is.na(x$relevance) && x$relevance >= min_relevance, + FUN.VALUE = base::logical(1) + )] + } + + cited_papers <- .hc_llm_unique_rag_papers(lapply(context_entries, function(x) x$paper)) + if (base::length(cited_papers) == 0) { + cited_papers <- .hc_llm_unique_rag_papers( + lapply(.hc_llm_rag_context_list(rag_result[["cited_papers"]]), .hc_llm_normalize_rag_paper) + ) + } + + out <- list( + query = .hc_llm_rag_scalar(rag_result[["query"]]), + answer = .hc_llm_rag_scalar(rag_result[["answer"]]), + context = context_entries, + cited_papers = cited_papers + ) + if (!base::is.null(rag_result[["raw_response_text"]])) { + out$raw_response_text <- .hc_llm_rag_scalar(rag_result[["raw_response_text"]]) + } + out +} + +.hc_llm_rag_context_list <- function(x) { + if (base::is.null(x)) { + return(list()) + } + if (base::is.data.frame(x)) { + return(lapply(base::seq_len(base::nrow(x)), function(i) { + row <- base::as.list(x[i, , drop = FALSE]) + lapply(row, function(v) { + if (base::length(v) == 1) v[[1]] else v + }) + })) + } + if (base::is.list(x) && base::length(x) > 0 && + any(base::names(x) %in% c("chunk", "section", "relevance", "paper", "apa_citation", "title", "doi"))) { + return(list(x)) + } + if (base::is.list(x)) { + return(x) + } + list(x) +} + +.hc_llm_normalize_rag_context_entry <- function(entry) { + if (!base::is.list(entry)) { + entry <- list(chunk = entry) + } + paper <- .hc_llm_normalize_rag_paper(entry[["paper"]]) + list( + chunk = .hc_llm_rag_scalar(entry[["chunk"]]), + section = .hc_llm_rag_scalar(entry[["section"]]), + relevance = .hc_llm_rag_numeric(entry[["relevance"]]), + paper = paper + ) +} + +.hc_llm_normalize_rag_paper <- function(paper) { + if (base::is.null(paper)) { + paper <- list() + } + if (base::is.data.frame(paper)) { + paper <- if (base::nrow(paper) > 0) { + base::as.list(paper[1, , drop = FALSE]) + } else { + list() + } + } + if (!base::is.list(paper)) { + paper <- list(apa_citation = paper) + } + list( + title = .hc_llm_rag_scalar(paper[["title"]]), + apa_citation = .hc_llm_rag_scalar(paper[["apa_citation"]]), + doi = .hc_llm_rag_scalar(paper[["doi"]]) + ) +} + +.hc_llm_rag_scalar <- function(x, default = "") { + if (base::is.null(x) || base::length(x) == 0) { + return(default) + } + if (base::is.data.frame(x)) { + if (base::nrow(x) == 0 || base::ncol(x) == 0) { + return(default) + } + x <- x[[1]][[1]] + } else if (base::is.list(x) && base::length(x) == 1 && !base::is.list(x[[1]])) { + x <- x[[1]] + } + val <- base::as.character(x[[1]]) + if (base::length(val) == 0 || base::is.na(val)) { + return(default) + } + base::trimws(val) +} + +.hc_llm_rag_numeric <- function(x) { + if (base::is.null(x) || base::length(x) == 0) { + return(NA_real_) + } + val <- suppressWarnings(base::as.numeric(x[[1]])) + if (base::length(val) == 0 || base::is.na(val)) { + return(NA_real_) + } + val[[1]] +} + +.hc_llm_unique_rag_papers <- function(papers) { + if (base::is.null(papers) || base::length(papers) == 0) { + return(list()) + } + papers <- lapply(papers, .hc_llm_normalize_rag_paper) + keys <- base::vapply( + papers, + function(p) base::paste(p$apa_citation, p$title, p$doi, sep = "\r"), + FUN.VALUE = base::character(1) + ) + has_value <- base::vapply( + papers, + function(p) base::any(base::nzchar(base::c(p$apa_citation, p$title, p$doi))), + FUN.VALUE = base::logical(1) + ) + keep <- !base::duplicated(keys) & has_value + papers[keep] +} + +.hc_llm_format_rag_context <- function(rag_result, max_chars = 12000) { + entries <- .hc_llm_rag_context_list(rag_result[["context"]]) + if (base::length(entries) == 0) { + return("") + } + + formatted_entries <- base::vapply( + base::seq_along(entries), + function(i) { + entry <- entries[[i]] + rel <- if (base::is.na(entry$relevance)) { + "NA" + } else { + base::format(base::round(entry$relevance, 3), nsmall = 3, trim = TRUE) + } + paper <- entry$paper + meta <- base::c( + base::paste0("[", i, "] relevance=", rel), + if (base::nzchar(entry$section)) base::paste0("section=", entry$section) else NULL, + if (base::nzchar(paper$title)) base::paste0("title=", paper$title) else NULL + ) + base::paste( + base::paste(meta, collapse = " | "), + if (base::nzchar(paper$apa_citation)) base::paste0("Citation: ", paper$apa_citation) else NULL, + if (base::nzchar(paper$doi)) base::paste0("DOI: ", paper$doi) else NULL, + base::paste0("Passage: ", entry$chunk), + sep = "\n" + ) + }, + FUN.VALUE = base::character(1) + ) + + citations <- .hc_llm_rag_citations_text(rag_result) + out <- base::paste( + "Use these passages only as supporting literature evidence. The gene list remains the primary evidence.", + base::paste(formatted_entries, collapse = "\n\n"), + if (base::nzchar(citations)) base::paste0("Unique cited papers: ", citations) else NULL, + sep = "\n\n" + ) + + if (!base::is.null(max_chars) && base::nchar(out, type = "chars") > max_chars) { + out <- base::paste0( + base::substr(out, 1L, max_chars), + "\n[DoRAG context truncated to ", max_chars, " characters.]" + ) + } + out +} + +.hc_llm_rag_context_count <- function(rag_result) { + base::length(.hc_llm_rag_context_list(rag_result[["context"]])) +} + +.hc_llm_rag_citations_text <- function(rag_result) { + papers <- .hc_llm_rag_context_list(rag_result[["cited_papers"]]) + if (base::length(papers) == 0) { + return("") + } + citations <- base::vapply( + papers, + function(p) { + p <- .hc_llm_normalize_rag_paper(p) + out <- if (base::nzchar(p$apa_citation)) p$apa_citation else p$title + if (base::nzchar(p$doi)) { + out <- base::paste0(out, " DOI: ", p$doi) + } + out + }, + FUN.VALUE = base::character(1) + ) + citations <- base::unique(citations[base::nzchar(citations)]) + base::paste(citations, collapse = " | ") +} + .hc_llm_error_result <- function(gene_info, label, context_text, @@ -1126,6 +1939,7 @@ hc_module_function_vllm <- function(...) { total_gene_count, context_text, truncated, + rag_context_text = NULL, llm = "gemini") { trunc_note <- if (isTRUE(truncated)) { base::paste0( @@ -1138,6 +1952,8 @@ hc_module_function_vllm <- function(...) { } biological_context <- if (base::nzchar(context_text)) context_text else "none provided" + has_rag_context <- !base::is.null(rag_context_text) && + base::nzchar(base::as.character(rag_context_text[[1]])) prompt_instructions <- if (llm == "vllm") { base::paste( @@ -1156,6 +1972,7 @@ hc_module_function_vllm <- function(...) { "Prefer specific states such as interferon-high inflammatory state, ribosome-high proliferative state, phagolysosomal activated state, antigen-presenting inflammatory state, platelet-like metabolic state, erythroid-skewed progenitor-like state, or macrophage-like transition state when supported.", "For `key_regulators`, list 2 to 5 likely transcription factors or signaling regulators separated by ' / '.", "If regulator evidence is weak, provide the most plausible regulators briefly rather than repeating the process.", + if (has_rag_context) "Use the retrieved DoRAG passages as optional literature support, but prioritize the supplied genes and biological context. Do not claim that RAG is statistical enrichment." else NULL, "Example style only:", '{"general_processes":"interferon signaling / antiviral innate immunity / antigen presentation","contextual_state":"activated interferon-high inflammatory monocyte state","key_regulators":"STAT1 / IRF7 / IRF9 / NFKB1"}' ) @@ -1173,7 +1990,8 @@ hc_module_function_vllm <- function(...) { "Provide `general_processes` as a short noun phrase listing the main biological program.", "Provide `contextual_state` as a short phrase describing the specific monocyte or transcriptional state in this study context.", "Provide `key_regulators` as a short phrase naming likely driving transcription factors or signaling regulators.", - "If regulator evidence is weak, state the most plausible regulators briefly rather than repeating the biological process." + "If regulator evidence is weak, state the most plausible regulators briefly rather than repeating the biological process.", + if (has_rag_context) "Use the retrieved DoRAG passages as optional literature support, but prioritize the supplied genes and biological context. Do not claim that RAG is statistical enrichment." else NULL ) } @@ -1183,6 +2001,14 @@ hc_module_function_vllm <- function(...) { base::paste0("Label: ", label), base::paste0("Gene-count note: ", trunc_note), base::paste0("Genes:\n", base::paste(genes, collapse = ", ")), + if (has_rag_context) { + base::paste0( + "Retrieved DoRAG literature context:\n", + base::as.character(rag_context_text[[1]]) + ) + } else { + NULL + }, sep = "\n" ) } @@ -1372,6 +2198,21 @@ hc_module_function_vllm <- function(...) { gene_count_input = base::integer(0), gene_count_sent = base::integer(0), truncated = base::logical(0), + rag_used = base::logical(0), + rag_query = base::character(0), + rag_context_count = base::integer(0), + rag_citations = base::character(0), + rag_error_message = base::character(0), + primary_interpretation_level = base::character(0), + without_context_general_processes = base::character(0), + without_context_contextual_state = base::character(0), + without_context_key_regulators = base::character(0), + with_context_general_processes = base::character(0), + with_context_contextual_state = base::character(0), + with_context_key_regulators = base::character(0), + with_rag_general_processes = base::character(0), + with_rag_contextual_state = base::character(0), + with_rag_key_regulators = base::character(0), status = base::character(0), error_message = base::character(0), timestamp = base::character(0), @@ -1397,6 +2238,21 @@ hc_module_function_vllm <- function(...) { gene_count_input = base::integer(0), gene_count_sent = base::integer(0), truncated = base::logical(0), + rag_used = base::logical(0), + rag_query = base::character(0), + rag_context_count = base::integer(0), + rag_citations = base::character(0), + rag_error_message = base::character(0), + primary_interpretation_level = base::character(0), + without_context_general_processes = base::character(0), + without_context_contextual_state = base::character(0), + without_context_key_regulators = base::character(0), + with_context_general_processes = base::character(0), + with_context_contextual_state = base::character(0), + with_context_key_regulators = base::character(0), + with_rag_general_processes = base::character(0), + with_rag_contextual_state = base::character(0), + with_rag_key_regulators = base::character(0), status = base::character(0), error_message = base::character(0), timestamp = base::character(0), @@ -1421,6 +2277,21 @@ hc_module_function_vllm <- function(...) { gene_count_input = .hc_llm_result_scalar(res, c("gene_count_input"), default = NA_integer_, mode = "integer"), gene_count_sent = .hc_llm_result_scalar(res, c("gene_count_sent"), default = NA_integer_, mode = "integer"), truncated = .hc_llm_result_scalar(res, c("truncated"), default = FALSE, mode = "logical"), + rag_used = .hc_llm_result_rag_used(res), + rag_query = .hc_llm_result_scalar(res, c("rag_query")), + rag_context_count = .hc_llm_result_rag_context_count(res), + rag_citations = .hc_llm_result_rag_citations(res), + rag_error_message = .hc_llm_result_scalar(res, c("rag_error_message")), + primary_interpretation_level = .hc_llm_result_scalar(res, c("primary_interpretation_level")), + without_context_general_processes = .hc_llm_result_level_scalar(res, "without_context", "general_processes"), + without_context_contextual_state = .hc_llm_result_level_scalar(res, "without_context", "contextual_state"), + without_context_key_regulators = .hc_llm_result_level_scalar(res, "without_context", "key_regulators"), + with_context_general_processes = .hc_llm_result_level_scalar(res, "with_context", "general_processes"), + with_context_contextual_state = .hc_llm_result_level_scalar(res, "with_context", "contextual_state"), + with_context_key_regulators = .hc_llm_result_level_scalar(res, "with_context", "key_regulators"), + with_rag_general_processes = .hc_llm_result_level_scalar(res, "with_rag", "general_processes"), + with_rag_contextual_state = .hc_llm_result_level_scalar(res, "with_rag", "contextual_state"), + with_rag_key_regulators = .hc_llm_result_level_scalar(res, "with_rag", "key_regulators"), status = .hc_llm_result_scalar(res, c("status")), error_message = .hc_llm_result_scalar(res, c("error_message")), timestamp = .hc_llm_result_scalar(res, c("timestamp")), @@ -1543,6 +2414,43 @@ hc_module_function_vllm <- function(...) { val } +.hc_llm_result_level_scalar <- function(res, level, field, default = NA_character_) { + val <- .hc_llm_result_field(res, c("interpretation_levels", level, "response", field)) + val <- .hc_llm_clean_text(val) + if (base::is.null(val) || base::length(val) == 0) { + return(default) + } + val <- base::as.character(val[[1]]) + if (base::length(val) == 0 || base::is.na(val) || !base::nzchar(val)) { + return(default) + } + val +} + +.hc_llm_result_rag_used <- function(res) { + !base::is.null(.hc_llm_result_field(res, c("rag"))) +} + +.hc_llm_result_rag_context_count <- function(res) { + rag <- .hc_llm_result_field(res, c("rag")) + if (base::is.null(rag)) { + return(0L) + } + base::as.integer(.hc_llm_rag_context_count(rag)) +} + +.hc_llm_result_rag_citations <- function(res) { + rag <- .hc_llm_result_field(res, c("rag")) + if (base::is.null(rag)) { + return(NA_character_) + } + citations <- .hc_llm_rag_citations_text(rag) + if (!base::nzchar(citations)) { + return(NA_character_) + } + citations +} + .hc_llm_short_title_fallback <- function(term, max_chars = 64) { if (base::is.null(term) || base::length(term) == 0) { return("No interpretation available") @@ -1637,6 +2545,21 @@ hc_module_function_vllm <- function(...) { gene_count_input = .hc_llm_result_scalar(res, c("gene_count_input"), default = NA_integer_, mode = "integer"), gene_count_sent = .hc_llm_result_scalar(res, c("gene_count_sent"), default = NA_integer_, mode = "integer"), truncated = .hc_llm_result_scalar(res, c("truncated"), default = FALSE, mode = "logical"), + rag_used = .hc_llm_result_rag_used(res), + rag_query = .hc_llm_result_scalar(res, c("rag_query")), + rag_context_count = .hc_llm_result_rag_context_count(res), + rag_citations = .hc_llm_result_rag_citations(res), + rag_error_message = .hc_llm_result_scalar(res, c("rag_error_message")), + primary_interpretation_level = .hc_llm_result_scalar(res, c("primary_interpretation_level")), + without_context_general_processes = .hc_llm_result_level_scalar(res, "without_context", "general_processes"), + without_context_contextual_state = .hc_llm_result_level_scalar(res, "without_context", "contextual_state"), + without_context_key_regulators = .hc_llm_result_level_scalar(res, "without_context", "key_regulators"), + with_context_general_processes = .hc_llm_result_level_scalar(res, "with_context", "general_processes"), + with_context_contextual_state = .hc_llm_result_level_scalar(res, "with_context", "contextual_state"), + with_context_key_regulators = .hc_llm_result_level_scalar(res, "with_context", "key_regulators"), + with_rag_general_processes = .hc_llm_result_level_scalar(res, "with_rag", "general_processes"), + with_rag_contextual_state = .hc_llm_result_level_scalar(res, "with_rag", "contextual_state"), + with_rag_key_regulators = .hc_llm_result_level_scalar(res, "with_rag", "key_regulators"), status = .hc_llm_result_scalar(res, c("status")), error_message = .hc_llm_result_scalar(res, c("error_message")), prompt = .hc_llm_result_scalar(res, c("prompt")), diff --git a/R/hc_plot_cluster_heatmap.R b/R/hc_plot_cluster_heatmap.R index 200ae41..68f6163 100644 --- a/R/hc_plot_cluster_heatmap.R +++ b/R/hc_plot_cluster_heatmap.R @@ -34,12 +34,15 @@ #' @param gene_count_mode Controls how gene counts per module are shown. One of "legacy", "bar_and_text", "bar", "text", "none". #' Default is "text". #' @param module_label_preset Layout preset for module labels. One of "auto", "compact", "balanced", "presentation". -#' Presets adjust automatic sizing only; explicitly set sizing parameters still take precedence. +#' Presets adjust automatic sizing only; explicitly set sizing parameters still +#' define the requested style, with a final fit guard applied at draw time. #' Default is "balanced". #' @param module_label_color Text color for labels drawn in module boxes. #' @param module_label_fontsize Optional numeric fontsize for module box labels. If NULL, a data-driven size is chosen. #' @param module_label_pt_size Optional numeric scale for label glyph size inside module boxes. #' Uses `snpc` units in `anno_simple()`. If NULL, a data-driven size is chosen. +#' The drawn size is automatically capped to one common safe size per heatmap +#' so all module names fit inside their boxes without mixed label sizes. #' @param module_box_width_cm Optional numeric width (cm) for module boxes. If NULL, width adapts to label length and fontsize. #' @param gene_count_fontsize Optional numeric fontsize for textual gene counts. If NULL, #' it defaults to the module label fontsize to keep both visually consistent. @@ -159,7 +162,8 @@ module_significance_show_qvalue = FALSE, module_significance_width_cm = 1.6, module_significance_p_cutoffs = c(0.001, 0.01, 0.05), - module_significance_annotation_name = "sig") { + module_significance_annotation_name = "sig", + write_module_tables = TRUE) { plot_cluster_heatmap_new( col_order = col_order, row_order = row_order, @@ -207,7 +211,8 @@ module_significance_show_qvalue = module_significance_show_qvalue, module_significance_width_cm = module_significance_width_cm, module_significance_p_cutoffs = module_significance_p_cutoffs, - module_significance_annotation_name = module_significance_annotation_name + module_significance_annotation_name = module_significance_annotation_name, + write_module_tables = write_module_tables ) } @@ -257,8 +262,15 @@ plot_cluster_heatmap <- function(col_order = NULL, module_significance_show_qvalue = FALSE, module_significance_width_cm = 1.6, module_significance_p_cutoffs = c(0.001, 0.01, 0.05), - module_significance_annotation_name = "sig") { + module_significance_annotation_name = "sig", + write_module_tables = TRUE) { .hc_alias_warning("plot_cluster_heatmap") + if (missing(module_label_numbering) && + .hc_module_label_map_has_split_labels( + hcobject[["integrated_output"]][["cluster_calc"]][["module_label_map"]] + )) { + module_label_numbering <- "preserve_existing" + } args <- as.list(environment()) invisible(base::do.call( .hc_run_modern_bridge, @@ -266,6 +278,18 @@ plot_cluster_heatmap <- function(col_order = NULL, )) } +.hc_module_label_map_has_split_labels <- function(module_label_map) { + if (base::is.null(module_label_map) || base::length(module_label_map) == 0) { + return(FALSE) + } + labels <- base::as.character(module_label_map) + labels <- labels[!base::is.na(labels) & base::nzchar(labels)] + # A split child carries a numeric suffix (e.g. "M3.1"); repeated splits append + # further suffixes ("M3.1.2"). Anchor on the trailing "." so nested + # splits are still detected. + base::any(base::grepl("\\.[0-9]+$", labels)) +} + .hc_module_label_draw_width_cm <- function(module_box_width_cm, module_labels_display = NULL, user_set_module_box_width_cm = FALSE, @@ -303,6 +327,384 @@ plot_cluster_heatmap <- function(col_order = NULL, draw_width } +.hc_cluster_heatmap_cell_size_mm <- function(n_heat_rows, + n_heat_cols, + duplicate_condition_width_scale = 1, + module_box_width_cm_draw = 0, + overall_plot_scale = 1) { + cell_size_mm <- 5.4 + if (n_heat_rows > 20) { + cell_size_mm <- 4.9 + } + if (n_heat_rows > 30) { + cell_size_mm <- 4.3 + } + if (n_heat_rows > 45) { + cell_size_mm <- 3.8 + } + if (n_heat_cols > 10) { + cell_size_mm <- base::min(cell_size_mm, 4.6) + } + if (n_heat_cols <= 4 && n_heat_rows <= 24) { + min_body_w_mm <- if (n_heat_cols <= 3) 30 else 34 + min_body_h_mm <- if (n_heat_rows <= 12) 90 else 108 + max_cell_mm <- if (n_heat_rows <= 12) 10 else 8 + boosted_cell_mm <- base::max( + min_body_w_mm / base::max(1, n_heat_cols), + min_body_h_mm / base::max(1, n_heat_rows) + ) + cell_size_mm <- base::max(cell_size_mm, base::min(max_cell_mm, boosted_cell_mm)) + } + if (duplicate_condition_width_scale > 1) { + target_module_box_to_cell_ratio <- 0.68 + min_cell_mm_from_box_ratio <- (module_box_width_cm_draw * 10) / target_module_box_to_cell_ratio + min_cell_mm_from_box_ratio <- base::min(10, min_cell_mm_from_box_ratio) + cell_size_mm <- base::max(cell_size_mm, min_cell_mm_from_box_ratio) + } + cell_size_mm * overall_plot_scale +} + +.hc_module_label_text_width_cm <- function(labels, + fontsize_pt, + fontface = "bold") { + labels_chr <- base::as.character(labels) + labels_chr[base::is.na(labels_chr)] <- "" + fontsize_pt <- base::as.numeric(fontsize_pt) + if (base::length(fontsize_pt) == 1) { + fontsize_pt <- base::rep(fontsize_pt, base::length(labels_chr)) + } + base::mapply( + FUN = function(label, font_pt) { + if (!base::nzchar(label) || !base::is.finite(font_pt) || font_pt <= 0) { + return(0) + } + tryCatch( + grid::convertWidth( + grid::grobWidth(grid::textGrob(label, gp = grid::gpar(fontsize = font_pt, fontface = fontface))), + "cm", + valueOnly = TRUE + ), + error = function(e) { + base::nchar(label) * font_pt * 0.021 + } + ) + }, + labels_chr, + fontsize_pt, + SIMPLIFY = TRUE, + USE.NAMES = FALSE + ) +} + +.hc_module_label_fit_pt <- function(module_label_pt_size, + module_box_width_cm, + module_labels_display = NULL, + n_heat_rows = 1, + cell_size_mm = 5.4, + module_label_fontsize = NULL, + use_fontsize_request = FALSE, + fontface = "bold") { + labels_chr <- base::as.character(module_labels_display) + labels_chr[base::is.na(labels_chr)] <- "" + n_labels <- base::length(labels_chr) + if (n_labels == 0) { + return(list( + pt_size = base::numeric(0), + base_pt_size = base::numeric(0), + available_width_cm = NA_real_, + available_height_pt = NA_real_, + width_limited = FALSE, + height_limited = FALSE, + shrunk = FALSE + )) + } + + module_label_pt_size <- .hc_first_numeric_value(module_label_pt_size) + module_label_fontsize <- .hc_first_numeric_value(module_label_fontsize) + module_box_width_cm <- .hc_first_numeric_value(module_box_width_cm) + n_heat_rows <- base::max(1L, base::as.integer(n_heat_rows)) + cell_size_mm <- .hc_first_numeric_value(cell_size_mm) + if (!base::is.finite(module_label_pt_size) || module_label_pt_size <= 0 || + !base::is.finite(module_box_width_cm) || module_box_width_cm <= 0 || + !base::is.finite(cell_size_mm) || cell_size_mm <= 0) { + return(list( + pt_size = base::rep(NA_real_, n_labels), + base_pt_size = base::rep(NA_real_, n_labels), + available_width_cm = NA_real_, + available_height_pt = NA_real_, + width_limited = FALSE, + height_limited = FALSE, + shrunk = FALSE + )) + } + + annotation_height_cm <- (n_heat_rows * cell_size_mm) / 10 + snpc_cm <- base::min(module_box_width_cm, annotation_height_cm) + base_pt <- if (isTRUE(use_fontsize_request) && + base::is.finite(module_label_fontsize) && + module_label_fontsize > 0) { + module_label_fontsize + } else { + module_label_pt_size * snpc_cm * 72 / 2.54 + } + base_pt_vec <- base::rep(base_pt, n_labels) + width_fill <- 0.95 + height_fill <- 0.95 + available_width_cm <- base::max(0.02, module_box_width_cm * width_fill) + available_height_pt <- (cell_size_mm * 72 / 25.4) * height_fill + + # `anno_simple(pch = "M1")` draws text-like symbols a little wider than a + # bare textGrob on some devices, so keep a tiny render guard while targeting + # 95% of the actual module box width. + symbol_width_guard <- 1.32 + width_at_10pt <- .hc_module_label_text_width_cm(labels_chr, fontsize_pt = 10, fontface = fontface) * + symbol_width_guard + max_pt_by_width <- base::ifelse( + width_at_10pt > 0, + 10 * available_width_cm / width_at_10pt, + Inf + ) + preferred_min_pt <- 2.8 + absolute_min_pt <- 0.2 + finite_width_pt <- max_pt_by_width[base::is.finite(max_pt_by_width)] + width_limit_pt <- if (base::length(finite_width_pt) > 0) { + base::min(finite_width_pt) + } else { + Inf + } + + common_pt <- base::min(base_pt, width_limit_pt, available_height_pt, na.rm = TRUE) + if (!base::is.finite(common_pt) || common_pt <= 0) { + common_pt <- base_pt + } + common_pt <- base::max(absolute_min_pt, common_pt) + if (common_pt >= preferred_min_pt || + (width_limit_pt >= preferred_min_pt && available_height_pt >= preferred_min_pt && base_pt >= preferred_min_pt)) { + common_pt <- base::max(preferred_min_pt, common_pt) + } + + fitted_pt <- base::rep(common_pt, n_labels) + fitted_width_cm <- .hc_module_label_text_width_cm(labels_chr, fontsize_pt = fitted_pt, fontface = fontface) * + symbol_width_guard + for (i in base::seq_len(5L)) { + max_width_cm <- base::max(fitted_width_cm, na.rm = TRUE) + if (!base::is.finite(max_width_cm) || + max_width_cm <= available_width_cm || + common_pt <= absolute_min_pt) { + break + } + common_pt <- base::pmax( + absolute_min_pt, + common_pt * (available_width_cm / max_width_cm) * 0.995 + ) + fitted_pt <- base::rep(common_pt, n_labels) + fitted_width_cm <- .hc_module_label_text_width_cm(labels_chr, fontsize_pt = fitted_pt, fontface = fontface) * + symbol_width_guard + } + + list( + pt_size = fitted_pt, + base_pt_size = base_pt_vec, + fitted_width_cm = fitted_width_cm, + available_width_cm = available_width_cm, + available_height_pt = available_height_pt, + width_limited = base::any(max_pt_by_width < base_pt_vec - 1e-8, na.rm = TRUE), + height_limited = base::any(available_height_pt < base_pt_vec - 1e-8, na.rm = TRUE), + below_preferred_min = base::is.finite(common_pt) && common_pt < preferred_min_pt, + shrunk = base::any(fitted_pt < base_pt_vec - 1e-8, na.rm = TRUE) + ) +} + +.hc_module_label_effective_fontsize <- function(module_label_fit, + fallback_fontsize = NULL, + fallback_pt_size = NULL, + min_pt = 0.2) { + fit_pt <- tryCatch( + .hc_as_numeric_safely(module_label_fit$pt_size), + error = function(e) numeric(0) + ) + fit_pt <- fit_pt[base::is.finite(fit_pt) & fit_pt > min_pt] + if (base::length(fit_pt) > 0) { + return(base::min(fit_pt)) + } + fallback <- .hc_first_numeric_value(fallback_fontsize) + if (base::is.finite(fallback) && fallback > min_pt) { + return(fallback) + } + fallback <- .hc_first_numeric_value(fallback_pt_size) + if (base::is.finite(fallback) && fallback > min_pt) { + return(fallback) + } + 5 +} + +.hc_module_label_box_annotation <- function(values, + colors, + labels = NULL, + label_color = "white", + label_fontsize_pt = NULL, + fontface = "bold", + width_cm, + border_gp = grid::gpar(col = "black", lwd = 0.5), + which = "row", + width_fill = 0.95, + height_fill = 0.90, + font_size_width_fill = 0.78, + min_visible_font_pt = 5.5, + minimum_sizing_label = "M4.2") { + values_chr <- base::as.character(values) + n_values <- base::length(values_chr) + labels_chr <- if (base::is.null(labels)) { + base::rep("", n_values) + } else { + label_names <- base::names(labels) + labels_raw <- base::as.character(labels) + if (base::length(labels_raw) == n_values) { + labels_raw + } else if (!base::is.null(label_names) && + base::length(label_names) == base::length(labels_raw) && + base::any(values_chr %in% base::as.character(label_names))) { + mapped <- base::as.character(labels_raw[base::match(values_chr, base::as.character(label_names))]) + missing_mapped <- base::is.na(mapped) | !base::nzchar(mapped) + mapped[missing_mapped] <- values_chr[missing_mapped] + mapped + } else if (base::length(labels_raw) > 0) { + base::rep_len(labels_raw, n_values) + } else { + values_chr + } + } + labels_chr[base::is.na(labels_chr)] <- "" + if (base::length(labels_chr) != n_values && n_values > 0) { + labels_chr <- base::rep_len(labels_chr, n_values) + } + + color_names <- base::names(colors) + colors_chr <- base::as.character(colors) + if (!base::is.null(color_names)) { + base::names(colors_chr) <- color_names + } + if (base::is.null(base::names(colors_chr))) { + base::names(colors_chr) <- base::as.character(colors_chr) + } + width_cm <- .hc_first_numeric_value(width_cm) + if (!base::is.finite(width_cm) || width_cm <= 0) { + width_cm <- 0.8 + } + label_fontsize_pt <- .hc_as_numeric_safely(label_fontsize_pt) + label_fontsize_pt <- label_fontsize_pt[base::is.finite(label_fontsize_pt) & label_fontsize_pt > 0] + max_font_pt <- if (base::length(label_fontsize_pt) > 0) { + base::min(label_fontsize_pt) + } else { + Inf + } + border_col <- if (!base::is.null(border_gp$col)) border_gp$col else "black" + border_lwd <- if (!base::is.null(border_gp$lwd)) border_gp$lwd else 0.5 + + fills_chr <- base::as.character(base::unname(colors_chr[values_chr])) + missing_fill <- base::is.na(fills_chr) | !base::nzchar(fills_chr) + fills_chr[missing_fill] <- values_chr[missing_fill] + valid_fill <- base::vapply(fills_chr, function(cl) { + base::isTRUE(tryCatch({ + grDevices::col2rgb(cl) + TRUE + }, error = function(e) FALSE)) + }, FUN.VALUE = base::logical(1)) + fills_chr[!valid_fill] <- "#d9d9d9" + + draw_fun <- function(index, ...) { + n <- base::length(index) + if (n == 0) { + return(invisible(NULL)) + } + y <- (n - base::seq_len(n) + 0.5) / n + fill <- fills_chr[index] + for (i in base::seq_len(n)) { + grid::grid.rect( + x = grid::unit(0.5, "npc"), + y = grid::unit(y[[i]], "npc"), + width = grid::unit(1, "npc"), + height = grid::unit(1 / n, "npc"), + gp = grid::gpar(fill = fill[[i]], col = border_col, lwd = border_lwd) + ) + } + + draw_labels <- labels_chr[index] + has_label <- !base::is.na(draw_labels) & base::nzchar(draw_labels) + if (!base::any(has_label)) { + return(invisible(NULL)) + } + all_labels <- labels_chr[!base::is.na(labels_chr) & base::nzchar(labels_chr)] + minimum_sizing_label <- base::as.character(minimum_sizing_label) + minimum_sizing_label <- minimum_sizing_label[!base::is.na(minimum_sizing_label) & + base::nzchar(minimum_sizing_label)] + sizing_labels <- base::unique(base::c(all_labels, minimum_sizing_label)) + target_width_cm <- grid::convertWidth(grid::unit(width_fill, "npc"), "cm", valueOnly = TRUE) + row_height_pt <- grid::convertHeight(grid::unit(1 / n, "npc"), "pt", valueOnly = TRUE) + target_height_pt <- base::max(row_height_pt * height_fill, min_visible_font_pt) + max_font_from_box_width_pt <- grid::convertWidth( + grid::unit(font_size_width_fill, "npc"), + "pt", + valueOnly = TRUE + ) + width_at_10pt <- .hc_module_label_text_width_cm(sizing_labels, fontsize_pt = 10, fontface = fontface) + max_width_10pt <- base::max(width_at_10pt, na.rm = TRUE) + width_fit_pt <- if (base::is.finite(max_width_10pt) && max_width_10pt > 0) { + 10 * target_width_cm / max_width_10pt + } else { + Inf + } + font_pt <- base::min( + max_font_pt, + width_fit_pt, + target_height_pt, + max_font_from_box_width_pt, + na.rm = TRUE + ) + if (!base::is.finite(font_pt) || font_pt <= 0) { + font_pt <- base::min(width_fit_pt, target_height_pt, max_font_from_box_width_pt, na.rm = TRUE) + } + if (!base::is.finite(font_pt) || font_pt <= 0) { + font_pt <- 5 + } + grid::grid.text( + draw_labels[has_label], + x = grid::unit(base::rep(0.5, base::sum(has_label)), "npc"), + y = grid::unit(y[has_label], "npc"), + gp = grid::gpar(col = label_color, fontsize = font_pt, fontface = fontface) + ) + invisible(NULL) + } + + ComplexHeatmap::AnnotationFunction( + fun = draw_fun, + fun_name = "module_label_boxes", + which = match.arg(which, c("column", "row")), + width = grid::unit(width_cm, "cm"), + n = n_values, + data_scale = c(0.5, 1.5), + var_import = list( + labels_chr = labels_chr, + fills_chr = fills_chr, + max_font_pt = max_font_pt, + label_color = label_color, + fontface = fontface, + width_fill = width_fill, + height_fill = height_fill, + font_size_width_fill = font_size_width_fill, + min_visible_font_pt = min_visible_font_pt, + minimum_sizing_label = minimum_sizing_label, + border_col = border_col, + border_lwd = border_lwd, + # `AnnotationFunction(var_import = ...)` re-homes `draw_fun` into an + # isolated environment whose parent is not the hcocena namespace, so the + # internal label-sizing helper must be imported explicitly or it cannot be + # found at draw time. + .hc_module_label_text_width_cm = .hc_module_label_text_width_cm + ) + ) +} + .hc_heatmap_screen_fit_scale <- function(total_width_mm, total_height_mm, device_size_mm = NULL, @@ -390,7 +792,8 @@ plot_cluster_heatmap_new <- function(col_order = NULL, module_significance_show_qvalue = FALSE, module_significance_width_cm = 1.6, module_significance_p_cutoffs = c(0.001, 0.01, 0.05), - module_significance_annotation_name = "sig") { + module_significance_annotation_name = "sig", + write_module_tables = TRUE) { pdf_width_is_default <- isTRUE(all.equal(as.numeric(pdf_width), 50)) pdf_height_is_default <- isTRUE(all.equal(as.numeric(pdf_height), 30)) @@ -551,6 +954,11 @@ plot_cluster_heatmap_new <- function(col_order = NULL, !base::nzchar(module_significance_annotation_name)) { stop("`module_significance_annotation_name` must be a non-empty string.") } + if (!base::is.logical(write_module_tables) || + base::length(write_module_tables) != 1 || + base::is.na(write_module_tables)) { + stop("`write_module_tables` must be TRUE or FALSE.") + } overall_plot_scale <- base::max(0.5, base::min(3, overall_plot_scale)) if (base::is.null(gfc_colors)) { gfc_colors <- .hc_default_gfc_colors() @@ -1337,6 +1745,12 @@ plot_cluster_heatmap_new <- function(col_order = NULL, base::max(base::nchar(module_labels_display), na.rm = TRUE) } + n_heat_rows <- base::nrow(mat_heatmap) + n_heat_cols <- base::ncol(mat_heatmap) + duplicate_condition_width_scale <- .hc_gfc_duplicate_condition_width_scale( + hcobject, + base::colnames(mat_heatmap) + ) n_rows <- base::nrow(c_df) preset_scale_font <- switch(module_label_preset, compact = 0.82, @@ -1371,7 +1785,7 @@ plot_cluster_heatmap_new <- function(col_order = NULL, base_width <- (max_chars * 0.09) + 0.15 } module_box_width_cm <- base::max( - 0.62, + 0.80, base::min( 4.8, base_width * preset_scale_width @@ -1385,25 +1799,33 @@ plot_cluster_heatmap_new <- function(col_order = NULL, module_sig_integrated = module_sig_integrated, max_sig_stars = max_sig_stars ) + cell_size_mm <- .hc_cluster_heatmap_cell_size_mm( + n_heat_rows = n_heat_rows, + n_heat_cols = n_heat_cols, + duplicate_condition_width_scale = duplicate_condition_width_scale, + module_box_width_cm_draw = module_box_width_cm_draw, + overall_plot_scale = overall_plot_scale + ) if (!user_set_module_label_pt_size) { base_pt <- if (n_rows <= 10) { - 0.45 + 0.95 } else if (n_rows <= 15) { - 0.30 + 0.90 } else if (n_rows <= 25) { - 0.26 + 0.82 } else if (n_rows <= 40) { - 0.22 + 0.72 } else { - 0.18 + 0.62 } - # No penalty for label length; box width adapts to character count instead + # Request a large default glyph; the box-fit guard chooses the largest + # common rendered size that still fits every module label. char_penalty <- 1 module_label_pt_size <- base::max( 0.14, base::min( - 0.55, + 1.05, (base_pt * preset_scale_pt) / char_penalty ) ) @@ -1416,10 +1838,33 @@ plot_cluster_heatmap_new <- function(col_order = NULL, base::min(0.55, module_label_pt_size * sig_pt_scale) ) } + module_label_fit <- .hc_module_label_fit_pt( + module_label_pt_size = module_label_pt_size_draw, + module_box_width_cm = module_box_width_cm_draw, + module_labels_display = module_labels_display, + n_heat_rows = n_heat_rows, + cell_size_mm = cell_size_mm, + module_label_fontsize = module_label_fontsize, + use_fontsize_request = user_set_module_label_fontsize && !user_set_module_label_pt_size, + fontface = "bold" + ) + module_label_pt_size_unit <- if (base::length(module_labels_display) > 0 && + base::length(module_label_fit$pt_size) == base::length(module_labels_display)) { + grid::unit(module_label_fit$pt_size, "pt") + } else { + grid::unit(module_label_pt_size_draw, "snpc") + } + module_label_fontsize_draw <- .hc_module_label_effective_fontsize( + module_label_fit = module_label_fit, + fallback_fontsize = module_label_fontsize, + fallback_pt_size = module_label_pt_size_draw + ) if (!user_set_gene_count_fontsize) { - # Keep numeric gene-count text aligned with module label size by default. - gene_count_fontsize <- module_label_fontsize + gene_count_fontsize <- base::max( + 5, + base::min(10, module_label_fontsize_draw * 0.72) + ) } if (!user_set_gene_count_pt_size) { @@ -1458,14 +1903,15 @@ plot_cluster_heatmap_new <- function(col_order = NULL, } else { module_significance_width_cm } - module_box_anno <- ComplexHeatmap::anno_simple( - row_order, - col = cluster_colors, - pch = module_labels_display, - pt_gp = grid::gpar(col = module_label_color, fontsize = module_label_fontsize, fontface = "bold"), - pt_size = grid::unit(module_label_pt_size_draw, "snpc"), - simple_anno_size = grid::unit(module_box_width_cm_draw, "cm"), - gp = module_box_border_gp, + module_box_anno <- .hc_module_label_box_annotation( + values = row_order, + colors = cluster_colors, + labels = module_labels_display, + label_color = module_label_color, + label_fontsize_pt = module_label_fit$base_pt_size, + fontface = "bold", + width_cm = module_box_width_cm_draw, + border_gp = module_box_border_gp, which = "row" ) @@ -1688,19 +2134,37 @@ plot_cluster_heatmap_new <- function(col_order = NULL, module_sig_width_cm_scaled + (0.8 * anno_scale) row_annotation_total_width_cm <- base_row_width_cm + enrichment_total_width_cm + module_label_fit_use <- .hc_module_label_fit_pt( + module_label_pt_size = module_label_pt_size_draw, + module_box_width_cm = module_box_width_cm_use, + module_labels_display = module_labels_display, + n_heat_rows = n_heat_rows, + cell_size_mm = cell_size_mm * anno_scale, + module_label_fontsize = module_label_fontsize, + use_fontsize_request = user_set_module_label_fontsize && !user_set_module_label_pt_size, + fontface = "bold" + ) + module_label_pt_size_unit_use <- if (base::length(module_labels_display) > 0 && + base::length(module_label_fit_use$pt_size) == base::length(module_labels_display)) { + grid::unit(module_label_fit_use$pt_size, "pt") + } else { + grid::unit(module_label_pt_size_draw, "snpc") + } + module_label_fontsize_draw_use <- .hc_module_label_effective_fontsize( + module_label_fit = module_label_fit_use, + fallback_fontsize = module_label_fontsize, + fallback_pt_size = module_label_pt_size_draw + ) - module_box_anno <- ComplexHeatmap::anno_simple( - row_order, - col = cluster_colors, - pch = module_labels_display, - pt_gp = grid::gpar( - col = module_label_color, - fontsize = module_label_fontsize, - fontface = "bold" - ), - pt_size = grid::unit(module_label_pt_size_draw, "snpc"), - simple_anno_size = grid::unit(module_box_width_cm_use, "cm"), - gp = module_box_border_gp, + module_box_anno <- .hc_module_label_box_annotation( + values = row_order, + colors = cluster_colors, + labels = module_labels_display, + label_color = module_label_color, + label_fontsize_pt = module_label_fit_use$base_pt_size, + fontface = "bold", + width_cm = module_box_width_cm_use, + border_gp = module_box_border_gp, which = "row" ) @@ -1967,10 +2431,6 @@ plot_cluster_heatmap_new <- function(col_order = NULL, column_labels_display <- .hc_gfc_display_col_labels(hcobject, base::colnames(mat_heatmap)) all_conditions <- .hc_gfc_display_count_labels(hcobject, base::colnames(mat_heatmap)) - duplicate_condition_width_scale <- .hc_gfc_duplicate_condition_width_scale( - hcobject, - base::colnames(mat_heatmap) - ) column_gap_k_enabled <- !base::is.numeric(k) || base::length(k) == 0 || base::all(k <= 0) column_gap_enabled <- (isTRUE(smart_column_gaps) || !base::is.null(column_gap_by)) && isTRUE(column_gap_k_enabled) @@ -1993,41 +2453,13 @@ plot_cluster_heatmap_new <- function(col_order = NULL, # --- 7. Plotting and Output --- # Keep module-expression tiles square-like for publication consistency. - n_heat_rows <- base::nrow(mat_heatmap) - n_heat_cols <- base::ncol(mat_heatmap) - cell_size_mm <- 5 - cell_size_mm <- 5.4 - if (n_heat_rows > 20) { - cell_size_mm <- 4.9 - } - if (n_heat_rows > 30) { - cell_size_mm <- 4.3 - } - if (n_heat_rows > 45) { - cell_size_mm <- 3.8 - } - if (n_heat_cols > 10) { - cell_size_mm <- base::min(cell_size_mm, 4.6) - } - if (n_heat_cols <= 4 && n_heat_rows <= 24) { - min_body_w_mm <- if (n_heat_cols <= 3) 30 else 34 - min_body_h_mm <- if (n_heat_rows <= 12) 90 else 108 - max_cell_mm <- if (n_heat_rows <= 12) 10 else 8 - boosted_cell_mm <- base::max( - min_body_w_mm / base::max(1, n_heat_cols), - min_body_h_mm / base::max(1, n_heat_rows) - ) - cell_size_mm <- base::max(cell_size_mm, base::min(max_cell_mm, boosted_cell_mm)) - } - if (duplicate_condition_width_scale > 1) { - # Keep module boxes visually narrower than one heatmap column when - # duplicated layer conditions add extra prefixed columns. - target_module_box_to_cell_ratio <- 0.68 - min_cell_mm_from_box_ratio <- (module_box_width_cm_draw * 10) / target_module_box_to_cell_ratio - min_cell_mm_from_box_ratio <- base::min(10, min_cell_mm_from_box_ratio) - cell_size_mm <- base::max(cell_size_mm, min_cell_mm_from_box_ratio) - } - cell_size_mm <- cell_size_mm * overall_plot_scale + cell_size_mm <- .hc_cluster_heatmap_cell_size_mm( + n_heat_rows = n_heat_rows, + n_heat_cols = n_heat_cols, + duplicate_condition_width_scale = duplicate_condition_width_scale, + module_box_width_cm_draw = module_box_width_cm_draw, + overall_plot_scale = overall_plot_scale + ) hm_width <- grid::unit((n_heat_cols * cell_size_mm) + column_gap_spec$total_gap_mm, "mm") hm_height <- grid::unit(n_heat_rows * cell_size_mm, "mm") max_col_chars <- if (base::length(column_labels_display) > 0) { @@ -2071,7 +2503,8 @@ plot_cluster_heatmap_new <- function(col_order = NULL, padding_obj = NULL, heatmap_legend_side = "right", annotation_legend_side = "right", - context = "cluster heatmap") { + context = "cluster heatmap", + capture_current_device = FALSE) { draw_once <- function(obj, lgd, show_ann_legend = TRUE, @@ -2094,18 +2527,47 @@ plot_cluster_heatmap_new <- function(col_order = NULL, } base::do.call(ComplexHeatmap::draw, args) } + run_attempt <- function(expr_fun) { + if (!isTRUE(capture_current_device)) { + return(list(result = expr_fun(), grob = NULL)) + } + result <- NULL + grob <- grid::grid.grabExpr( + { + result <- expr_fun() + }, + wrap = TRUE + ) + list(result = result, grob = grob) + } + finish_attempt <- function(attempt) { + if (isTRUE(capture_current_device) && !base::is.null(attempt$grob)) { + grid::grid.newpage() + grid::grid.draw(attempt$grob) + } + attempt$result + } - first <- try(draw_once(ht_obj, legend_obj, show_ann_legend = TRUE), silent = TRUE) + first <- try( + run_attempt(function() draw_once(ht_obj, legend_obj, show_ann_legend = TRUE)), + silent = TRUE + ) if (!inherits(first, "try-error")) { - return(first) + return(finish_attempt(first)) } - try(grid::grid.newpage(), silent = TRUE) + if (!isTRUE(capture_current_device)) { + try(grid::grid.newpage(), silent = TRUE) + } second <- try( - draw_once( - deep_clone(ht_obj), - deep_clone(legend_obj), - show_ann_legend = TRUE + run_attempt( + function() { + draw_once( + deep_clone(ht_obj), + deep_clone(legend_obj), + show_ann_legend = TRUE + ) + } ), silent = TRUE ) @@ -2116,15 +2578,21 @@ plot_cluster_heatmap_new <- function(col_order = NULL, ".", call. = FALSE ) - return(second) + return(finish_attempt(second)) } - try(grid::grid.newpage(), silent = TRUE) + if (!isTRUE(capture_current_device)) { + try(grid::grid.newpage(), silent = TRUE) + } third <- try( - draw_once( - deep_clone(ht_obj), - NULL, - show_ann_legend = FALSE + run_attempt( + function() { + draw_once( + deep_clone(ht_obj), + NULL, + show_ann_legend = FALSE + ) + } ), silent = TRUE ) @@ -2135,17 +2603,23 @@ plot_cluster_heatmap_new <- function(col_order = NULL, ": annotation legends were disabled after viewport issues.", call. = FALSE ) - return(third) + return(finish_attempt(third)) } - try(grid::grid.newpage(), silent = TRUE) + if (!isTRUE(capture_current_device)) { + try(grid::grid.newpage(), silent = TRUE) + } fourth <- try( - draw_once( - deep_clone(ht_obj), - NULL, - show_ann_legend = FALSE, - show_heat_legend = TRUE, - use_padding = FALSE + run_attempt( + function() { + draw_once( + deep_clone(ht_obj), + NULL, + show_ann_legend = FALSE, + show_heat_legend = TRUE, + use_padding = FALSE + ) + } ), silent = TRUE ) @@ -2156,17 +2630,23 @@ plot_cluster_heatmap_new <- function(col_order = NULL, ": legends/padding were simplified after viewport issues.", call. = FALSE ) - return(fourth) + return(finish_attempt(fourth)) } - try(grid::grid.newpage(), silent = TRUE) + if (!isTRUE(capture_current_device)) { + try(grid::grid.newpage(), silent = TRUE) + } fifth <- try( - draw_once( - deep_clone(ht_obj), - NULL, - show_ann_legend = FALSE, - show_heat_legend = FALSE, - use_padding = FALSE + run_attempt( + function() { + draw_once( + deep_clone(ht_obj), + NULL, + show_ann_legend = FALSE, + show_heat_legend = FALSE, + use_padding = FALSE + ) + } ), silent = TRUE ) @@ -2177,7 +2657,7 @@ plot_cluster_heatmap_new <- function(col_order = NULL, ": all legends were disabled after viewport issues.", call. = FALSE ) - return(fifth) + return(finish_attempt(fifth)) } stop(fifth) @@ -2424,7 +2904,8 @@ plot_cluster_heatmap_new <- function(col_order = NULL, padding_obj = draw_padding_screen, heatmap_legend_side = heatmap_legend_side_mode, annotation_legend_side = annotation_legend_side_mode, - context = "cluster heatmap current device" + context = "cluster heatmap current device", + capture_current_device = TRUE ), error = function(e) { warning( @@ -2487,32 +2968,110 @@ plot_cluster_heatmap_new <- function(col_order = NULL, stringsAsFactors = FALSE ) } - module_gene_list_file <- base::paste0( - hcobject[["working_directory"]][["dir_output"]], - hcobject[["global_settings"]][["save_folder"]], - "/Module_Gene_List.xlsx" - ) - tryCatch( + if (base::isTRUE(write_module_tables)) { + module_gene_list_file <- base::paste0( + hcobject[["working_directory"]][["dir_output"]], + hcobject[["global_settings"]][["save_folder"]], + "/Module_Gene_List.xlsx" + ) + tryCatch( + { + openxlsx::write.xlsx( + x = list(module_gene_list = module_gene_list_tbl), + file = module_gene_list_file, + overwrite = TRUE + ) + }, + error = function(e) { + warning( + "Could not write Module_Gene_List.xlsx: ", + base::conditionMessage(e) + ) + } + ) + } + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_map"), module_label_map) + .hc_set_bridge_hcobject_slot(c("satellite_outputs", "module_gene_list"), module_gene_list_tbl) + + # Mean-GFC-per-module-and-group table: the numeric values behind the heatmap + # cells (rows = modules, columns = sample groups). `mat_heatmap` is keyed by + # cluster colour, so relabel rows with the displayed module labels and keep the + # colour as a reference column. + module_gfc_means_tbl <- tryCatch( { - openxlsx::write.xlsx( - x = list(module_gene_list = module_gene_list_tbl), - file = module_gene_list_file, - overwrite = TRUE + ordered_colors <- base::as.character(row_order) + ordered_colors <- ordered_colors[ordered_colors %in% base::rownames(mat_heatmap)] + ordered_colors <- base::c( + ordered_colors, + base::setdiff(base::rownames(mat_heatmap), ordered_colors) + ) + gfc_means_mat <- mat_heatmap[ordered_colors, , drop = FALSE] + row_labels <- base::as.character(module_export_map[ordered_colors]) + missing_lab <- base::is.na(row_labels) | !base::nzchar(row_labels) + row_labels[missing_lab] <- ordered_colors[missing_lab] + base::data.frame( + module = row_labels, + cluster_color = ordered_colors, + base::as.data.frame(gfc_means_mat, check.names = FALSE, stringsAsFactors = FALSE), + check.names = FALSE, + stringsAsFactors = FALSE, + row.names = NULL ) }, error = function(e) { - warning( - "Could not write Module_Gene_List.xlsx: ", - base::conditionMessage(e) + base::warning( + "Could not build module GFC means table: ", + base::conditionMessage(e), + call. = FALSE ) + NULL } ) - .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_map"), module_label_map) - .hc_set_bridge_hcobject_slot(c("satellite_outputs", "module_gene_list"), module_gene_list_tbl) + if (!base::is.null(module_gfc_means_tbl)) { + if (base::isTRUE(write_module_tables)) { + module_gfc_means_file <- base::paste0( + hcobject[["working_directory"]][["dir_output"]], + hcobject[["global_settings"]][["save_folder"]], + "/Module_GFC_Means.xlsx" + ) + tryCatch( + openxlsx::write.xlsx( + x = base::list(module_gfc_means = module_gfc_means_tbl), + file = module_gfc_means_file, + overwrite = TRUE + ), + error = function(e) { + base::warning( + "Could not write Module_GFC_Means.xlsx: ", + base::conditionMessage(e) + ) + } + ) + } + .hc_set_bridge_hcobject_slot(c("satellite_outputs", "module_gfc_means"), module_gfc_means_tbl) + } + + module_label_fit_pt <- .hc_as_numeric_safely(module_label_fit$pt_size) + module_label_fit_base_pt <- .hc_as_numeric_safely(module_label_fit$base_pt_size) + module_label_fit_pt <- module_label_fit_pt[base::is.finite(module_label_fit_pt)] + module_label_fit_base_pt <- module_label_fit_base_pt[base::is.finite(module_label_fit_base_pt)] + module_label_fit_min <- if (base::length(module_label_fit_pt) > 0) base::min(module_label_fit_pt) else NA_real_ + module_label_fit_max <- if (base::length(module_label_fit_pt) > 0) base::max(module_label_fit_pt) else NA_real_ + module_label_fit_requested_max <- if (base::length(module_label_fit_base_pt) > 0) { + base::max(module_label_fit_base_pt) + } else { + NA_real_ + } .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_mode"), module_label_mode) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_numbering"), module_label_numbering) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_fontsize"), module_label_fontsize) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_pt_size"), module_label_pt_size) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_pt_size_effective_pt_min"), module_label_fit_min) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_pt_size_effective_pt_max"), module_label_fit_max) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_pt_size_requested_pt_max"), module_label_fit_requested_max) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_auto_fit_shrunk"), isTRUE(module_label_fit$shrunk)) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_auto_fit_width_limited"), isTRUE(module_label_fit$width_limited)) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_label_auto_fit_height_limited"), isTRUE(module_label_fit$height_limited)) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "module_box_width_cm"), module_box_width_cm) module_box_to_cell_ratio <- (.hc_first_numeric_value(module_box_width_cm) * 10) / .hc_first_numeric_value(cell_size_mm) diff --git a/R/hc_plot_cluster_heatmap_view.R b/R/hc_plot_cluster_heatmap_view.R new file mode 100644 index 0000000..321de2d --- /dev/null +++ b/R/hc_plot_cluster_heatmap_view.R @@ -0,0 +1,404 @@ +#' Plot a fixed-scale view of selected modules from the module heatmap +#' +#' Creates a module heatmap view for selected modules in the exact order +#' supplied by `modules`. The module labels and module-color annotations are +#' preserved, and the GFC expression color scale is inherited from the full +#' heatmap cache instead of being recomputed from the subset. +#' +#' @param hc A `HCoCenaExperiment`. +#' @param modules Character/numeric vector of modules to show. Accepts current +#' module labels such as `"M2"` or `"M2.1"`, module colors, or numeric module +#' indices in the current heatmap order. +#' @param file_name Optional export file name for the view. Use `FALSE` to skip +#' file export. +#' @param col_order Optional condition order. If `NULL`, the stored full +#' heatmap column order is reused when available. +#' @param cluster_columns Logical. Whether to cluster columns in the view. +#' @param cluster_rows Logical. Whether to cluster rows in the view. Defaults to +#' `FALSE` so the `modules` vector controls the displayed order. +#' @param gfc_colors Optional GFC palette override. If `NULL`, the stored full +#' heatmap palette is reused. +#' @param gfc_scale_limits Optional GFC scale override. If `NULL`, the stored +#' full heatmap scale is reused; if unavailable, a full-cache/global fallback +#' is used. +#' @param preserve_full_heatmap Logical. If `TRUE`, restore the pre-existing +#' full heatmap cache in `hc` after writing the view. +#' @param store_view Logical. If `TRUE`, store lightweight metadata for the +#' view under `hc@integration@cluster$heatmap_views`. +#' @param view_name Optional name used when `store_view = TRUE`. +#' @param ... Additional plotting arguments forwarded to +#' [hc_plot_cluster_heatmap()]. +#' +#' @return Updated `HCoCenaExperiment`. +#' @export +hc_plot_cluster_heatmap_view <- function(hc, + modules, + file_name = "Heatmap_modules_view.pdf", + col_order = NULL, + cluster_columns = FALSE, + cluster_rows = FALSE, + gfc_colors = NULL, + gfc_scale_limits = NULL, + preserve_full_heatmap = TRUE, + store_view = TRUE, + view_name = NULL, + ...) { + if (!inherits(hc, "HCoCenaExperiment")) { + stop("`hc` must be a `HCoCenaExperiment`.") + } + if (missing(modules) || base::is.null(modules) || base::length(modules) == 0) { + stop("`modules` must contain at least one module label, color, or index.") + } + if (!base::is.logical(cluster_columns) || base::length(cluster_columns) != 1 || + base::is.na(cluster_columns)) { + stop("`cluster_columns` must be TRUE or FALSE.") + } + if (!base::is.logical(cluster_rows) || base::length(cluster_rows) != 1 || + base::is.na(cluster_rows)) { + stop("`cluster_rows` must be TRUE or FALSE.") + } + if (!base::is.logical(preserve_full_heatmap) || + base::length(preserve_full_heatmap) != 1 || + base::is.na(preserve_full_heatmap)) { + stop("`preserve_full_heatmap` must be TRUE or FALSE.") + } + if (!base::is.logical(store_view) || + base::length(store_view) != 1 || + base::is.na(store_view)) { + stop("`store_view` must be TRUE or FALSE.") + } + + dot_args <- list(...) + if ("row_order" %in% base::names(dot_args)) { + stop("Do not pass `row_order`; use `modules` to control the view order.") + } + if ("write_module_tables" %in% base::names(dot_args)) { + stop("`write_module_tables` is managed by `hc_plot_cluster_heatmap_view()`.") + } + + cluster_calc <- .hc_heatmap_view_cluster_calc(hc) + resolved_modules <- .hc_heatmap_view_resolve_modules( + modules = modules, + cluster_calc = cluster_calc + ) + if (base::is.null(col_order)) { + heatmap_info <- tryCatch(.hc_heatmap_cache_info(cluster_calc), error = function(e) NULL) + if (!base::is.null(heatmap_info) && + !base::is.null(heatmap_info$col_order) && + base::length(heatmap_info$col_order) > 0) { + col_order <- base::as.character(heatmap_info$col_order) + } + } + + resolved_gfc_colors <- .hc_heatmap_view_resolve_gfc_colors( + cluster_calc = cluster_calc, + gfc_colors = gfc_colors + ) + resolved_gfc_scale_limits <- .hc_heatmap_view_resolve_gfc_scale_limits( + hc = hc, + cluster_calc = cluster_calc, + gfc_scale_limits = gfc_scale_limits + ) + + plot_args <- list( + hc = hc, + file_name = file_name, + row_order = resolved_modules$target_colors, + col_order = col_order, + cluster_rows = cluster_rows, + cluster_columns = cluster_columns, + gfc_colors = resolved_gfc_colors, + gfc_scale_limits = resolved_gfc_scale_limits, + write_module_tables = FALSE + ) + if (!("module_label_numbering" %in% base::names(dot_args))) { + plot_args[["module_label_numbering"]] <- "preserve_existing" + } + plot_args <- base::c(plot_args, dot_args) + + plotted_hc <- base::do.call(hc_plot_cluster_heatmap, plot_args) + view_state <- .hc_heatmap_view_capture_state( + hc = plotted_hc, + requested_modules = modules, + resolved_modules = resolved_modules, + file_name = file_name, + gfc_colors = resolved_gfc_colors, + gfc_scale_limits = resolved_gfc_scale_limits + ) + + out <- if (base::isTRUE(preserve_full_heatmap)) { + .hc_heatmap_view_restore_full_state( + plotted_hc = plotted_hc, + original_hc = hc + ) + } else { + plotted_hc + } + + if (base::isTRUE(store_view)) { + out <- .hc_heatmap_view_store_state( + hc = out, + view_name = view_name, + view_state = view_state + ) + } + + methods::validObject(out) + out +} + +.hc_heatmap_view_cluster_calc <- function(hc) { + cluster_calc <- tryCatch(hc@integration@cluster, error = function(e) NULL) + if (base::is.null(cluster_calc)) { + return(list()) + } + base::as.list(cluster_calc) +} + +.hc_heatmap_view_available_modules <- function(cluster_calc) { + cluster_info <- tryCatch(cluster_calc[["cluster_information"]], error = function(e) NULL) + if (base::is.data.frame(cluster_info) && + "color" %in% base::colnames(cluster_info) && + base::nrow(cluster_info) > 0) { + included <- base::rep(TRUE, base::nrow(cluster_info)) + if ("cluster_included" %in% base::colnames(cluster_info)) { + included <- base::as.character(cluster_info$cluster_included) == "yes" + included[base::is.na(included)] <- FALSE + } + colors <- base::as.character(cluster_info$color[included]) + colors <- colors[!base::is.na(colors) & base::nzchar(colors)] + if (base::length(colors) > 0) { + return(base::unique(colors)) + } + } + + heatmap_info <- tryCatch(.hc_heatmap_cache_info(cluster_calc), error = function(e) NULL) + if (!base::is.null(heatmap_info) && + !base::is.null(heatmap_info$matrix) && + !base::is.null(base::rownames(heatmap_info$matrix))) { + row_ids <- base::as.character(base::rownames(heatmap_info$matrix)) + row_ids <- row_ids[!base::is.na(row_ids) & base::nzchar(row_ids)] + if (base::length(row_ids) > 0) { + return(base::unique(row_ids)) + } + } + + module_label_map <- tryCatch(cluster_calc[["module_label_map"]], error = function(e) NULL) + module_names <- base::names(module_label_map) + module_names <- base::as.character(module_names) + module_names <- module_names[!base::is.na(module_names) & base::nzchar(module_names)] + base::unique(module_names) +} + +.hc_heatmap_view_resolve_modules <- function(modules, cluster_calc) { + available_colors <- .hc_heatmap_view_available_modules(cluster_calc) + if (base::length(available_colors) == 0) { + stop( + "No module identifiers are available. Run `hc_plot_cluster_heatmap()` ", + "or check that `hc@integration@cluster` contains module information." + ) + } + + module_prefix <- tryCatch(cluster_calc[["module_prefix"]], error = function(e) NULL) + module_prefix <- base::as.character(module_prefix) + if (base::length(module_prefix) != 1 || + base::is.na(module_prefix) || + !base::nzchar(module_prefix)) { + module_prefix <- "M" + } + module_label_map <- .hc_normalize_module_label_map_for_split( + module_label_map = tryCatch(cluster_calc[["module_label_map"]], error = function(e) NULL), + available_colors = available_colors, + module_prefix = module_prefix + ) + module_order <- .hc_split_module_order( + cluster_calc = cluster_calc, + available_colors = available_colors + ) + resolved <- .hc_resolve_modules_for_split( + modules = modules, + available_colors = available_colors, + module_label_map = module_label_map, + module_order = module_order + ) + + unresolved_tbl <- resolved$resolution_table[ + base::as.character(resolved$resolution_table$status) != "ok", , + drop = FALSE + ] + if (base::nrow(unresolved_tbl) > 0) { + unresolved_inputs <- base::unique(base::as.character(unresolved_tbl$input)) + available_labels <- base::unique(base::as.character(module_label_map[module_order])) + available_labels <- available_labels[!base::is.na(available_labels) & base::nzchar(available_labels)] + if (base::length(available_labels) > 0) { + available_preview <- base::paste(utils::head(available_labels, 12L), collapse = ", ") + if (base::length(available_labels) > 12L) { + available_preview <- base::paste0(available_preview, ", ...") + } + stop( + "Could not resolve the following module identifier(s): ", + base::paste(unresolved_inputs, collapse = ", "), + "\nUse exact current module labels, module colors, or numeric indices.", + "\nAvailable module labels: ", available_preview + ) + } + stop( + "Could not resolve the following module identifier(s): ", + base::paste(unresolved_inputs, collapse = ", "), + "\nUse exact current module labels, module colors, or numeric indices." + ) + } + if (base::length(resolved$target_colors) == 0) { + stop( + "Could not match any requested module in `modules`.\n", + "Use module labels, module colors, or module indices." + ) + } + + resolved +} + +.hc_heatmap_view_normalize_gfc_scale_limits <- function(x) { + if (base::is.null(x)) { + return(NULL) + } + x <- .hc_as_numeric_safely(x) + if (base::length(x) == 1) { + if (!base::is.finite(x) || x <= 0) { + stop("`gfc_scale_limits` as single value must be finite and > 0.") + } + return(c(-base::abs(x), base::abs(x))) + } + if (base::length(x) != 2 || base::any(!base::is.finite(x))) { + stop("`gfc_scale_limits` must be NULL, one positive number, or a numeric vector of length 2.") + } + x <- base::sort(x) + if (x[[1]] == x[[2]]) { + stop("`gfc_scale_limits` must have different min/max values.") + } + x +} + +.hc_heatmap_view_resolve_gfc_scale_limits <- function(hc, + cluster_calc, + gfc_scale_limits = NULL) { + resolved <- .hc_heatmap_view_normalize_gfc_scale_limits(gfc_scale_limits) + if (!base::is.null(resolved)) { + return(resolved) + } + + stored <- tryCatch( + .hc_heatmap_view_normalize_gfc_scale_limits(cluster_calc[["gfc_scale_limits"]]), + error = function(e) NULL + ) + if (!base::is.null(stored)) { + return(stored) + } + + heatmap_info <- tryCatch(.hc_heatmap_cache_info(cluster_calc), error = function(e) NULL) + mat <- if (!base::is.null(heatmap_info)) heatmap_info$matrix else NULL + if (!base::is.null(mat)) { + mat_num <- base::suppressWarnings(base::as.numeric(mat)) + mat_num <- mat_num[base::is.finite(mat_num)] + if (base::length(mat_num) > 0) { + lim <- base::max(base::abs(mat_num), na.rm = TRUE) + if (base::is.finite(lim) && lim > 0) { + return(c(-base::abs(lim), base::abs(lim))) + } + } + } + + bridge <- tryCatch(.hc_as_bridge_object_for_cluster_plot(hc), error = function(e) NULL) + fallback_lim <- tryCatch( + .hc_first_numeric_value(bridge[["global_settings"]][["range_GFC"]]), + error = function(e) NA_real_ + ) + if (base::length(fallback_lim) != 1 || + !base::is.finite(fallback_lim) || + fallback_lim <= 0) { + fallback_lim <- 2 + } + c(-base::abs(fallback_lim), base::abs(fallback_lim)) +} + +.hc_heatmap_view_resolve_gfc_colors <- function(cluster_calc, gfc_colors = NULL) { + colors <- gfc_colors + if (base::is.null(colors)) { + colors <- tryCatch(base::as.character(cluster_calc[["gfc_colors"]]), error = function(e) NULL) + } + if (base::is.null(colors) || + base::length(colors) < 2 || + !base::is.character(colors) || + base::any(base::is.na(colors)) || + base::any(colors == "")) { + colors <- .hc_default_gfc_colors() + } + base::as.character(colors) +} + +.hc_heatmap_view_capture_state <- function(hc, + requested_modules, + resolved_modules, + file_name, + gfc_colors, + gfc_scale_limits) { + cluster_calc <- .hc_heatmap_view_cluster_calc(hc) + keep_names <- c( + "heatmap_output_files", + "heatmap_matrix", + "heatmap_row_order", + "heatmap_column_order", + "heatmap_column_labels_display", + "heatmap_cluster", + "heatmap_cluster_raw", + "module_label_map", + "module_label_mode", + "module_label_numbering", + "gfc_colors", + "gfc_scale_limits", + "overall_plot_scale" + ) + view_state <- cluster_calc[base::intersect(keep_names, base::names(cluster_calc))] + view_state[["requested_modules"]] <- requested_modules + view_state[["resolved_modules"]] <- resolved_modules$target_colors + view_state[["resolved_labels"]] <- resolved_modules$resolved_labels + view_state[["resolution_table"]] <- resolved_modules$resolution_table + view_state[["file_name"]] <- file_name + view_state[["gfc_colors"]] <- gfc_colors + view_state[["gfc_scale_limits"]] <- gfc_scale_limits + view_state[["created_at"]] <- base::Sys.time() + view_state +} + +.hc_heatmap_view_restore_full_state <- function(plotted_hc, original_hc) { + plotted_hc@integration@cluster <- original_hc@integration@cluster + plotted_hc@satellite <- original_hc@satellite + plotted_hc +} + +.hc_heatmap_view_store_state <- function(hc, view_name = NULL, view_state) { + if (base::is.null(view_name)) { + view_name <- base::paste0("view_", base::format(base::Sys.time(), "%Y%m%d_%H%M%S")) + } + if (!base::is.character(view_name) || + base::length(view_name) != 1 || + base::is.na(view_name) || + !base::nzchar(base::trimws(view_name))) { + stop("`view_name` must be NULL or a non-empty character scalar.") + } + view_name <- base::trimws(view_name) + + cluster_calc <- .hc_heatmap_view_cluster_calc(hc) + views <- tryCatch(cluster_calc[["heatmap_views"]], error = function(e) NULL) + if (base::is.null(views)) { + views <- list() + } else { + views <- base::as.list(views) + } + views[[view_name]] <- view_state + cluster_calc[["heatmap_views"]] <- views + + hc@integration@cluster <- base::do.call(S4Vectors::SimpleList, cluster_calc) + methods::validObject(hc) + hc +} diff --git a/R/hc_plot_module_function_gemini.R b/R/hc_plot_module_function_gemini.R index 41c7848..a98152f 100644 --- a/R/hc_plot_module_function_gemini.R +++ b/R/hc_plot_module_function_gemini.R @@ -10,8 +10,14 @@ #' `"llm_module_function"`. #' @param modules Optional character vector to subset modules. #' @param fields Character vector selecting which LLM fields to plot. Supported -#' values are `"general_processes"`, `"contextual_state"`, and -#' `"key_regulators"`. Defaults to all three. +#' values include `"general_processes"`, `"contextual_state"`, +#' `"key_regulators"`, and comparison fields from +#' `compare_interpretation_levels = TRUE`, such as +#' `"without_context_contextual_state"`, `"with_context_contextual_state"`, +#' and `"with_rag_contextual_state"`. The aliases +#' `"contextual_state_without_context"`, `"contextual_state_with_context"`, +#' and `"contextual_state_rag"` are also accepted. Defaults to the three +#' top-level fields. #' @param max_chars Maximum number of characters shown per term. Default is #' `90`. #' @param text_size Numeric text size passed to `ggplot2::geom_text()`. @@ -158,9 +164,10 @@ hc_plot_module_function_llm <- function(hc, if (!is.numeric(dpi) || length(dpi) != 1 || is.na(dpi) || dpi <= 0) { stop("`dpi` must be a single positive number.") } + field_map <- .hc_llm_plot_field_map() fields <- unique(as.character(fields)) fields <- fields[!is.na(fields) & fields != ""] - valid_fields <- c("general_processes", "contextual_state", "key_regulators") + valid_fields <- names(field_map$source) if (length(fields) == 0 || !all(fields %in% valid_fields)) { stop("`fields` must contain one or more of: ", paste(valid_fields, collapse = ", "), ".") } @@ -195,7 +202,8 @@ hc_plot_module_function_llm <- function(hc, summary_tbl$module <- as.character(summary_tbl$module) summary_tbl$module_color <- as.character(summary_tbl$module_color) summary_tbl$module_color[is.na(summary_tbl$module_color) | summary_tbl$module_color == ""] <- "grey70" - for (nm in valid_fields) { + needed_summary_fields <- unique(unname(field_map$source[fields])) + for (nm in needed_summary_fields) { if (!nm %in% colnames(summary_tbl)) { summary_tbl[[nm]] <- NA_character_ } @@ -205,26 +213,21 @@ hc_plot_module_function_llm <- function(hc, summary_tbl$label_color <- vapply(summary_tbl$module_color, .hc_llm_contrast_text_color, FUN.VALUE = character(1)) summary_tbl$module_factor <- factor(summary_tbl$module, levels = rev(summary_tbl$module)) - title_map <- c( - general_processes = "General processes", - contextual_state = "Contextual state", - key_regulators = "Key regulators" - ) - out <- lapply(fields, function(field_nm) { field_tbl <- summary_tbl + source_nm <- unname(field_map$source[[field_nm]]) field_tbl$term_plot <- vapply( ifelse( - is.na(field_tbl[[field_nm]]) | field_tbl[[field_nm]] == "", + is.na(field_tbl[[source_nm]]) | field_tbl[[source_nm]] == "", "No interpretation available.", - field_tbl[[field_nm]] + field_tbl[[source_nm]] ), .hc_llm_prepare_display_title, FUN.VALUE = character(1), max_chars = as.integer(max_chars[[1]]) ) - plot_title <- paste0(title, ": ", title_map[[field_nm]]) + plot_title <- paste0(title, ": ", field_map$title[[field_nm]]) if (isTRUE(with_heatmap) && !is.null(heatmap_info) && isTRUE(heatmap_info$draw_supported)) { combined_grob <- .hc_llm_capture_combined_heatmap_grob( @@ -335,6 +338,56 @@ hc_plot_module_function_gemini <- function(...) { hc_plot_module_function_llm(...) } +.hc_llm_plot_field_map <- function() { + source <- c( + general_processes = "general_processes", + contextual_state = "contextual_state", + key_regulators = "key_regulators", + without_context_general_processes = "without_context_general_processes", + without_context_contextual_state = "without_context_contextual_state", + without_context_key_regulators = "without_context_key_regulators", + with_context_general_processes = "with_context_general_processes", + with_context_contextual_state = "with_context_contextual_state", + with_context_key_regulators = "with_context_key_regulators", + with_rag_general_processes = "with_rag_general_processes", + with_rag_contextual_state = "with_rag_contextual_state", + with_rag_key_regulators = "with_rag_key_regulators", + general_processes_without_context = "without_context_general_processes", + contextual_state_without_context = "without_context_contextual_state", + key_regulators_without_context = "without_context_key_regulators", + general_processes_with_context = "with_context_general_processes", + contextual_state_with_context = "with_context_contextual_state", + key_regulators_with_context = "with_context_key_regulators", + general_processes_rag = "with_rag_general_processes", + contextual_state_rag = "with_rag_contextual_state", + key_regulators_rag = "with_rag_key_regulators" + ) + title <- c( + general_processes = "General processes", + contextual_state = "Contextual state", + key_regulators = "Key regulators", + without_context_general_processes = "General processes without context", + without_context_contextual_state = "Contextual state without context", + without_context_key_regulators = "Key regulators without context", + with_context_general_processes = "General processes with context", + with_context_contextual_state = "Contextual state with context", + with_context_key_regulators = "Key regulators with context", + with_rag_general_processes = "General processes with RAG", + with_rag_contextual_state = "Contextual state with RAG", + with_rag_key_regulators = "Key regulators with RAG", + general_processes_without_context = "General processes without context", + contextual_state_without_context = "Contextual state without context", + key_regulators_without_context = "Key regulators without context", + general_processes_with_context = "General processes with context", + contextual_state_with_context = "Contextual state with context", + key_regulators_with_context = "Key regulators with context", + general_processes_rag = "General processes with RAG", + contextual_state_rag = "Contextual state with RAG", + key_regulators_rag = "Key regulators with RAG" + ) + list(source = source, title = title) +} + .hc_llm_plot_output_dir <- function(hc) { paths <- tryCatch(.hc_row_to_list(hc@config@paths), error = function(e) list()) global_cfg <- tryCatch(.hc_row_to_list(hc@config@global), error = function(e) list()) @@ -679,7 +732,7 @@ print.hc_llm_heatmap_plot <- function(x, ...) { } module_label_fontsize <- 4.2 module_box_width_cm <- base::max(0.56, base::min(4.8, (max_label_chars * 0.09) + 0.15)) - module_label_pt_size <- 0.42 + module_label_pt_size <- 0.90 cell_size_mm <- .hc_llm_default_heatmap_cell_size_mm( n_heat_rows = n_heat_rows, n_heat_cols = n_heat_cols, @@ -953,15 +1006,38 @@ print.hc_llm_heatmap_plot <- function(x, ...) { grid_width = grid::unit(3.2 * overall_plot_scale, "mm"), legend_height = grid::unit(legend_height_mm, "mm") ) + module_label_fit <- .hc_module_label_fit_pt( + module_label_pt_size = module_pt_size, + module_box_width_cm = module_box_width_cm, + module_labels_display = module_by_row, + n_heat_rows = nrow(mat_use), + cell_size_mm = cell_mm, + module_label_fontsize = label_fontsize, + use_fontsize_request = !base::is.null(module_label_fontsize) && + base::is.null(module_label_pt_size), + fontface = "bold" + ) + module_label_pt_size_unit <- if (base::length(module_by_row) > 0 && + base::length(module_label_fit$pt_size) == base::length(module_by_row)) { + grid::unit(module_label_fit$pt_size, "pt") + } else { + grid::unit(module_pt_size, "snpc") + } + module_label_fontsize_draw <- .hc_module_label_effective_fontsize( + module_label_fit = module_label_fit, + fallback_fontsize = label_fontsize, + fallback_pt_size = module_pt_size + ) - module_box_anno <- ComplexHeatmap::anno_simple( - module_by_row, - col = module_col_map, - pch = module_by_row, - pt_gp = grid::gpar(col = "white", fontsize = label_fontsize, fontface = "bold"), - pt_size = grid::unit(module_pt_size, "snpc"), - simple_anno_size = grid::unit(module_box_width_cm, "cm"), - gp = module_box_border_gp, + module_box_anno <- .hc_module_label_box_annotation( + values = module_by_row, + colors = module_col_map, + labels = module_by_row, + label_color = "white", + label_fontsize_pt = module_label_fit$base_pt_size, + fontface = "bold", + width_cm = module_box_width_cm, + border_gp = module_box_border_gp, which = "row" ) diff --git a/R/hc_split_modules.R b/R/hc_split_modules.R index ea18fee..9e2b591 100644 --- a/R/hc_split_modules.R +++ b/R/hc_split_modules.R @@ -9,7 +9,8 @@ #' Accepted values in `modules`: #' - module labels from `module_label_map` (e.g. `"M3"`), #' - module colors (e.g. `"#FFD700"`), -#' - numeric module indices (based on current included-module order). +#' - numeric module indices (based on the current heatmap row order when +#' available, otherwise the current included-module order). #' #' The function stores an undo snapshot in #' `hcobject$satellite_outputs$module_split_history`. @@ -20,7 +21,10 @@ #' `"cluster_fast_greedy"`, `"cluster_infomap"`, `"cluster_walktrap"`, #' `"cluster_label_prop"` or `"auto"`. #' @param no_of_iterations Number of Leiden iterations (used only for Leiden). -#' @param resolution Leiden resolution (used only for Leiden). +#' @param resolution Leiden resolution (used only for Leiden). Use either one +#' positive value for all selected modules or one positive value per selected +#' module in the same order as `modules`. Named vectors may use module labels +#' or module colors. #' @param resolution_grid Optional numeric vector of candidate resolutions to #' test before splitting. For each candidate, hCoCena reports how many #' submodules would be retained after size filtering. @@ -60,8 +64,11 @@ if (!is.numeric(no_of_iterations) || length(no_of_iterations) != 1 || !is.finite(no_of_iterations) || no_of_iterations < 1) { stop("`no_of_iterations` must be a positive numeric scalar.") } - if (!is.numeric(resolution) || length(resolution) != 1 || !is.finite(resolution) || resolution <= 0) { - stop("`resolution` must be a positive numeric scalar.") + if (!is.numeric(resolution) || + length(resolution) == 0 || + any(!is.finite(resolution)) || + any(resolution <= 0)) { + stop("`resolution` must contain positive finite numeric value(s).") } if (!is.null(resolution_grid) && (!is.numeric(resolution_grid) || length(resolution_grid) == 0 || any(!is.finite(resolution_grid)) || any(resolution_grid <= 0))) { @@ -70,6 +77,9 @@ if (!is.logical(resolution_test_only) || length(resolution_test_only) != 1 || is.na(resolution_test_only)) { stop("`resolution_test_only` must be TRUE or FALSE.") } + if (isTRUE(resolution_test_only) && is.null(resolution_grid)) { + stop("`resolution_test_only = TRUE` requires `resolution_grid`; no split was applied.") + } if (!is.character(partition_type) || length(partition_type) != 1 || is.na(partition_type)) { stop("`partition_type` must be a single character string.") } @@ -135,21 +145,23 @@ "M" } - module_label_map <- cluster_calc[["module_label_map"]] - if (is.null(module_label_map) || length(module_label_map) == 0) { - module_label_map <- stats::setNames( - paste0(module_prefix, seq_len(nrow(incl_info))), - as.character(incl_info$color) - ) - } else { - module_label_map <- as.character(module_label_map) - names(module_label_map) <- as.character(names(cluster_calc[["module_label_map"]])) - } + raw_module_label_map <- cluster_calc[["module_label_map"]] + module_label_map <- raw_module_label_map + module_label_map <- .hc_normalize_module_label_map_for_split( + module_label_map = module_label_map, + available_colors = as.character(incl_info$color), + module_prefix = module_prefix + ) + module_order <- .hc_split_module_order( + cluster_calc = cluster_calc, + available_colors = as.character(incl_info$color) + ) resolved <- .hc_resolve_modules_for_split( modules = modules, available_colors = as.character(incl_info$color), - module_label_map = module_label_map + module_label_map = module_label_map, + module_order = module_order ) unresolved_tbl <- resolved$resolution_table[ as.character(resolved$resolution_table$status) != "ok", , @@ -157,7 +169,7 @@ ] if (nrow(unresolved_tbl) > 0) { unresolved_inputs <- unique(as.character(unresolved_tbl$input)) - available_labels <- unique(as.character(module_label_map)) + available_labels <- unique(as.character(module_label_map[module_order])) if (length(available_labels) > 0) { available_preview <- paste(utils::head(available_labels, 12L), collapse = ", ") if (length(available_labels) > 12L) { @@ -183,6 +195,21 @@ "Use module labels (e.g. M3), module colors, or module indices." ) } + resolution_by_module <- .hc_split_resolution_by_module( + resolution = resolution, + target_colors = target_colors, + resolved_labels = resolved$resolved_labels, + resolution_table = resolved$resolution_table + ) + resolved$resolution_table$resolution <- NA_real_ + resolution_match <- match( + as.character(resolved$resolution_table$resolved_color), + names(resolution_by_module) + ) + has_resolution_match <- !is.na(resolution_match) + resolved$resolution_table$resolution[has_resolution_match] <- as.numeric( + resolution_by_module[resolution_match[has_resolution_match]] + ) if (isTRUE(verbose)) { message( @@ -190,6 +217,7 @@ paste(resolved$resolved_labels, collapse = ", "), " (", length(target_colors), " module(s))." ) + .hc_display_object(resolved$resolution_table, row.names = FALSE) } all_graph_nodes <- as.character(igraph::V(merged_net)$name) @@ -223,6 +251,8 @@ modules = as.character(modules), cluster_algo = cluster_algo, no_of_iterations = as.integer(round(no_of_iterations)), + resolution = resolution, + resolution_by_module = resolution_by_module, resolution_grid = resolution_grid, partition_type = partition_type, seed = as.integer(round(seed)), @@ -258,7 +288,11 @@ } before_cluster_info <- cluster_info - before_module_label_map <- module_label_map + before_module_label_map <- if (!is.null(raw_module_label_map)) { + raw_module_label_map + } else { + module_label_map + } new_rows <- vector("list", nrow(cluster_info)) row_ptr <- 0L @@ -290,11 +324,13 @@ genes <- unique(genes) if (length(genes) < 3) { + this_resolution <- as.numeric(resolution_by_module[[this_color]]) row_ptr <- row_ptr + 1L new_rows[[row_ptr]] <- row_i split_summary_rows[[length(split_summary_rows) + 1L]] <- data.frame( parent_color = this_color, parent_label = parent_label, + resolution = this_resolution, parent_genes = length(genes), n_submodules_raw = 1L, n_submodules = 1L, @@ -307,11 +343,12 @@ } subgraph <- igraph::induced_subgraph(merged_net, vids = genes) + this_resolution <- as.numeric(resolution_by_module[[this_color]]) membership <- .hc_split_membership( graph_obj = subgraph, cluster_algo = cluster_algo, no_of_iterations = as.integer(round(no_of_iterations)), - resolution = resolution, + resolution = this_resolution, partition_type = partition_type, seed = as.integer(round(seed)) ) @@ -322,6 +359,7 @@ split_summary_rows[[length(split_summary_rows) + 1L]] <- data.frame( parent_color = this_color, parent_label = parent_label, + resolution = this_resolution, parent_genes = length(genes), n_submodules_raw = 1L, n_submodules = 1L, @@ -355,6 +393,7 @@ split_summary_rows[[length(split_summary_rows) + 1L]] <- data.frame( parent_color = this_color, parent_label = parent_label, + resolution = this_resolution, parent_genes = length(genes), n_submodules_raw = n_sub_raw, n_submodules = 0L, @@ -370,17 +409,23 @@ split_summary_rows[[length(split_summary_rows) + 1L]] <- data.frame( parent_color = this_color, parent_label = parent_label, + resolution = this_resolution, parent_genes = length(genes), n_submodules_raw = n_sub_raw, n_submodules = 1L, - removed_small_submodules = removed_small_submodules, - removed_small_genes = removed_small_genes, + removed_small_submodules = 0L, + removed_small_genes = 0L, status = if (isTRUE(drop_small_submodules)) "skipped_after_small_module_filter" else "skipped_single_submodule", stringsAsFactors = FALSE ) next } member_levels <- member_levels[kept_idx] + member_levels <- .hc_order_split_member_levels( + membership = membership, + member_levels = member_levels, + gfc_all = gfc_all + ) n_sub <- length(member_levels) total_removed_small_submodules <- total_removed_small_submodules + removed_small_submodules total_removed_small_genes <- total_removed_small_genes + removed_small_genes @@ -414,6 +459,7 @@ split_summary_rows[[length(split_summary_rows) + 1L]] <- data.frame( parent_color = this_color, parent_label = parent_label, + resolution = this_resolution, parent_genes = length(genes), n_submodules_raw = n_sub_raw, n_submodules = n_sub, @@ -482,6 +528,7 @@ cluster_algo = cluster_algo, no_of_iterations = as.integer(round(no_of_iterations)), resolution = resolution, + resolution_by_module = resolution_by_module, partition_type = partition_type, seed = as.integer(round(seed)), drop_small_submodules = drop_small_submodules, @@ -777,14 +824,96 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { )) } -.hc_resolve_modules_for_split <- function(modules, available_colors, module_label_map) { +.hc_normalize_module_label_map_for_split <- function(module_label_map, + available_colors, + module_prefix = "M") { available_colors <- unique(as.character(available_colors)) + if (length(available_colors) == 0) { + return(character()) + } + + if (is.null(module_label_map) || length(module_label_map) == 0) { + return(stats::setNames( + paste0(module_prefix, seq_along(available_colors)), + available_colors + )) + } + + map_names <- names(module_label_map) + module_label_map <- as.character(module_label_map) + if (!is.null(map_names) && length(map_names) == length(module_label_map)) { + names(module_label_map) <- as.character(map_names) + } + + missing_before <- setdiff(available_colors, names(module_label_map)) + if (length(missing_before) > 0 && !is.null(names(module_label_map))) { + inverse_map <- stats::setNames(names(module_label_map), as.character(module_label_map)) + if (all(available_colors %in% names(inverse_map))) { + module_label_map <- inverse_map + } + } + + module_label_map <- module_label_map[names(module_label_map) %in% available_colors] + missing_map <- setdiff(available_colors, names(module_label_map)) + if (length(missing_map) > 0) { + start_idx <- length(module_label_map) + 1L + module_label_map <- c( + module_label_map, + stats::setNames( + paste0(module_prefix, seq.int(start_idx, length.out = length(missing_map))), + missing_map + ) + ) + } + + module_label_map[available_colors] +} + +.hc_split_module_order <- function(cluster_calc, available_colors) { + available_colors <- unique(as.character(available_colors)) + if (length(available_colors) == 0) { + return(character()) + } + + row_order <- tryCatch(cluster_calc[["heatmap_row_order"]], error = function(e) NULL) + row_order <- as.character(unlist(row_order, use.names = FALSE)) + row_order <- row_order[!is.na(row_order) & nzchar(row_order) & row_order %in% available_colors] + if (length(row_order) == 0) { + heatmap_info <- tryCatch(.hc_heatmap_cache_info(cluster_calc), error = function(e) NULL) + row_order <- if (!is.null(heatmap_info)) { + as.character(heatmap_info$row_order) + } else { + character() + } + row_order <- row_order[!is.na(row_order) & nzchar(row_order) & row_order %in% available_colors] + } + + if (length(row_order) == 0) { + return(available_colors) + } + unique(c(row_order, setdiff(available_colors, row_order))) +} + +.hc_resolve_modules_for_split <- function(modules, + available_colors, + module_label_map, + module_order = NULL) { + available_colors <- unique(as.character(available_colors)) + module_order <- unique(as.character(module_order)) + module_order <- module_order[!is.na(module_order) & nzchar(module_order) & module_order %in% available_colors] + if (length(module_order) == 0) { + module_order <- available_colors + } else { + module_order <- c(module_order, setdiff(available_colors, module_order)) + } label_to_color <- stats::setNames(names(module_label_map), as.character(module_label_map)) resolved_colors <- character(0) resolution_rows <- list() input_vec <- modules - if (is.numeric(input_vec)) { + if (is.list(input_vec)) { + input_vec <- as.list(input_vec) + } else if (is.numeric(input_vec)) { input_vec <- as.list(input_vec) } else { input_vec <- as.list(as.character(input_vec)) @@ -797,9 +926,9 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { resolved <- NA_character_ if (is.numeric(x) && is.finite(x)) { - idx <- as.integer(round(x)) - if (idx >= 1 && idx <= length(available_colors)) { - resolved <- available_colors[[idx]] + idx <- if (abs(x - round(x)) < 1e-8) as.integer(round(x)) else NA_integer_ + if (!is.na(idx) && idx >= 1 && idx <= length(module_order)) { + resolved <- module_order[[idx]] } else { status <- "not_found" } @@ -809,6 +938,13 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { resolved <- x_chr } else if (x_chr %in% names(label_to_color)) { resolved <- as.character(label_to_color[[x_chr]]) + } else if (grepl("^[0-9]+$", x_chr)) { + idx <- suppressWarnings(as.integer(x_chr)) + if (!is.na(idx) && idx >= 1 && idx <= length(module_order)) { + resolved <- module_order[[idx]] + } else { + status <- "not_found" + } } else { status <- "not_found" } @@ -828,9 +964,15 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { } else { NA_character_ } + resolved_index <- if (!is.na(resolved)) { + match(resolved, module_order) + } else { + NA_integer_ + } resolution_rows[[length(resolution_rows) + 1L]] <- data.frame( input = x_chr, + resolved_index = as.integer(resolved_index), resolved_color = resolved, resolved_label = resolved_label, status = status, @@ -855,6 +997,130 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { ) } +.hc_split_resolution_by_module <- function(resolution, + target_colors, + resolved_labels, + resolution_table = NULL) { + if (!is.numeric(resolution) || + length(resolution) == 0 || + any(!is.finite(resolution)) || + any(resolution <= 0)) { + stop("`resolution` must contain positive finite numeric value(s).") + } + + target_colors <- unique(as.character(target_colors)) + resolved_labels <- as.character(resolved_labels) + if (length(target_colors) == 0) { + return(stats::setNames(numeric(0), character(0))) + } + if (length(resolved_labels) != length(target_colors)) { + resolved_labels <- target_colors + } + + resolution_names <- names(resolution) + resolution <- as.numeric(resolution) + if (!is.null(resolution_names) && length(resolution_names) == length(resolution)) { + names(resolution) <- as.character(resolution_names) + } + + resolution_names <- names(resolution) + resolution_names_trimmed <- if (is.null(resolution_names)) { + character(length(resolution)) + } else { + trimws(as.character(resolution_names)) + } + has_any_names <- any(!is.na(resolution_names_trimmed) & nzchar(resolution_names_trimmed)) + if (!has_any_names) { + if (length(resolution) == 1) { + return(stats::setNames(rep(as.numeric(resolution[[1]]), length(target_colors)), target_colors)) + } + if (length(resolution) != length(target_colors)) { + stop( + "`resolution` must be length 1 or have one value per resolved module (", + length(target_colors), "); got ", length(resolution), "." + ) + } + return(stats::setNames(as.numeric(resolution), target_colors)) + } + + resolution_names <- resolution_names_trimmed + if (length(resolution_names) != length(resolution) || + any(is.na(resolution_names) | !nzchar(resolution_names))) { + stop("All values in a named `resolution` vector must have non-empty names.") + } + + input_to_color <- character(0) + if (!is.null(resolution_table) && + is.data.frame(resolution_table) && + all(c("input", "resolved_color") %in% colnames(resolution_table))) { + ok <- !is.na(resolution_table$resolved_color) + if (any(ok)) { + input_to_color <- stats::setNames( + as.character(resolution_table$resolved_color[ok]), + as.character(resolution_table$input[ok]) + ) + } + } + + out <- stats::setNames(rep(NA_real_, length(target_colors)), target_colors) + unmatched <- character(0) + duplicated_targets <- character(0) + + for (i in seq_along(resolution)) { + nm <- resolution_names[[i]] + matched_color <- character(0) + + if (nm %in% target_colors) { + matched_color <- nm + } else { + label_hits <- target_colors[resolved_labels == nm] + input_hits <- input_to_color[names(input_to_color) == nm] + matched_color <- unique(c(label_hits, unname(input_hits))) + matched_color <- matched_color[matched_color %in% target_colors] + } + + if (length(matched_color) == 0) { + unmatched <- c(unmatched, nm) + next + } + if (length(matched_color) > 1) { + duplicated_targets <- c(duplicated_targets, nm) + next + } + + color <- matched_color[[1]] + if (!is.na(out[[color]]) && !isTRUE(all.equal(out[[color]], as.numeric(resolution[[i]])))) { + duplicated_targets <- c(duplicated_targets, nm) + next + } + out[[color]] <- as.numeric(resolution[[i]]) + } + + if (length(unmatched) > 0) { + stop( + "`resolution` names must match requested module labels, module colors, or input module identifiers. ", + "Unmatched: ", paste(unique(unmatched), collapse = ", "), "." + ) + } + if (length(duplicated_targets) > 0) { + stop( + "`resolution` names must map to one unique requested module. Ambiguous or conflicting names: ", + paste(unique(duplicated_targets), collapse = ", "), "." + ) + } + + missing_idx <- which(is.na(out)) + if (length(missing_idx) > 0) { + missing_desc <- paste0(resolved_labels[missing_idx], " (", target_colors[missing_idx], ")") + stop( + "`resolution` must provide one value for each resolved module. Missing: ", + paste(missing_desc, collapse = ", "), "." + ) + } + + out +} + .hc_split_membership <- function(graph_obj, cluster_algo = "cluster_leiden", no_of_iterations = 2L, @@ -962,6 +1228,49 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { out } +.hc_order_split_member_levels <- function(membership, + member_levels, + gfc_all) { + member_levels <- as.character(member_levels) + if (length(member_levels) <= 1 || is.null(gfc_all) || !is.data.frame(gfc_all) || !("Gene" %in% colnames(gfc_all))) { + return(member_levels) + } + + # Submodules are ordered by hierarchical clustering of their mean GFC + # profiles. Euclidean distance is invariant to column permutation, so the + # displayed heatmap column order has no bearing on this ordering. + value_idx <- .hc_gfc_value_col_idx(gfc_all) + if (length(value_idx) == 0) { + return(member_levels) + } + + row_means <- lapply(member_levels, function(lvl) { + genes_k <- names(membership)[membership == as.integer(lvl)] + genes_k <- unique(genes_k[!is.na(genes_k) & genes_k != ""]) + sub <- gfc_all[gfc_all$Gene %in% genes_k, value_idx, drop = FALSE] + if (nrow(sub) == 0) { + return(rep(NA_real_, length(value_idx))) + } + mat <- as.matrix(data.frame(lapply(sub, identity), check.names = FALSE)) + storage.mode(mat) <- "numeric" + colMeans(mat, na.rm = TRUE) + }) + profile_mat <- do.call(rbind, row_means) + rownames(profile_mat) <- member_levels + profile_mat[!is.finite(profile_mat)] <- 0 + + if (nrow(profile_mat) <= 1 || ncol(profile_mat) == 0) { + return(member_levels) + } + ordered <- tryCatch({ + rownames(profile_mat)[stats::hclust(stats::dist(profile_mat), method = "complete")$order] + }, error = function(e) NULL) + if (is.null(ordered) || length(ordered) == 0) { + return(member_levels) + } + unique(c(as.character(ordered), setdiff(member_levels, ordered))) +} + .hc_build_child_cluster_rows <- function(template_row, membership, member_levels, @@ -969,7 +1278,6 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { child_labels, gfc_all) { out <- template_row[rep(1, length(member_levels)), , drop = FALSE] - gfc_cols <- setdiff(colnames(gfc_all), "Gene") for (k in seq_along(member_levels)) { lvl <- member_levels[[k]] @@ -996,19 +1304,12 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { out$vertexsize[[k]] <- 3 } - if (length(gfc_cols) > 0) { - if ("conditions" %in% colnames(out)) { - out$conditions[[k]] <- paste0(gfc_cols, collapse = "#") - } - if ("grp_means" %in% colnames(out)) { - sub <- gfc_all[gfc_all$Gene %in% genes_k, gfc_cols, drop = FALSE] - means <- if (nrow(sub) > 0) { - colMeans(sub, na.rm = TRUE) - } else { - rep(NA_real_, length(gfc_cols)) - } - out$grp_means[[k]] <- paste0(round(means, 3), collapse = ",") - } + gfc_means <- .hc_gfc_colmeans_for_genes(gfc_all, genes = genes_k) + if ("conditions" %in% colnames(out)) { + out$conditions[[k]] <- paste0(names(gfc_means), collapse = "#") + } + if ("grp_means" %in% colnames(out)) { + out$grp_means[[k]] <- paste0(round(gfc_means, 3), collapse = ",") } } @@ -1077,6 +1378,10 @@ unsplit_modules <- function(which = c("last", "all"), verbose = TRUE) { .hc_set_bridge_hcobject_slot(c("integrated_output", "knowledge_network"), NULL) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "labelled_network"), NULL) .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "network_col_by_module"), NULL) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "heatmap_cluster"), NULL) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "heatmap_cluster_raw"), NULL) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "heatmap_matrix"), NULL) + .hc_set_bridge_hcobject_slot(c("integrated_output", "cluster_calc", "heatmap_row_order"), NULL) sat <- hcobject[["satellite_outputs"]] if (is.null(sat) || !is.list(sat)) { diff --git a/R/hc_transcription_factor_enrichment.R b/R/hc_transcription_factor_enrichment.R index b4cd9d1..83444fa 100644 --- a/R/hc_transcription_factor_enrichment.R +++ b/R/hc_transcription_factor_enrichment.R @@ -268,6 +268,7 @@ TF_overrep_module <- function(clusters = "all", topTF = 5, topTarget = 5) { page_labels = page_labels, width = 15, height = 8, + display = TRUE, draw_page_fun = function(page_idx, page_label) { module_idx <- ((page_idx - 1L) %/% 2L) + 1L if ((page_idx %% 2L) == 1L) { diff --git a/R/hc_upstream_inference.R b/R/hc_upstream_inference.R index e4ca2d1..1f50e83 100644 --- a/R/hc_upstream_inference.R +++ b/R/hc_upstream_inference.R @@ -2291,14 +2291,36 @@ upstream_inference <- function(resources = c("TF", "Pathway"), } else { base::max(0.22, base::min(0.90, 0.075 * label_fontsize)) } - module_box_anno <- ComplexHeatmap::anno_simple( - keep_clusters, - col = module_color_map, - pch = module_labels, - pt_gp = grid::gpar(col = "white", fontsize = label_fontsize, fontface = "bold"), - pt_size = grid::unit(module_label_pt_size, "snpc"), - simple_anno_size = grid::unit(module_box_width_cm, "cm"), - gp = grid::gpar(col = "black"), + module_label_fit <- .hc_module_label_fit_pt( + module_label_pt_size = module_label_pt_size, + module_box_width_cm = module_box_width_cm, + module_labels_display = module_labels, + n_heat_rows = n_hc_rows, + cell_size_mm = hc_cell_mm, + module_label_fontsize = label_fontsize, + use_fontsize_request = FALSE, + fontface = "bold" + ) + module_label_pt_size_unit <- if (base::length(module_labels) > 0 && + base::length(module_label_fit$pt_size) == base::length(module_labels)) { + grid::unit(module_label_fit$pt_size, "pt") + } else { + grid::unit(module_label_pt_size, "snpc") + } + module_label_fontsize_draw <- .hc_module_label_effective_fontsize( + module_label_fit = module_label_fit, + fallback_fontsize = label_fontsize, + fallback_pt_size = module_label_pt_size + ) + module_box_anno <- .hc_module_label_box_annotation( + values = keep_clusters, + colors = module_color_map, + labels = module_labels, + label_color = "white", + label_fontsize_pt = module_label_fit$base_pt_size, + fontface = "bold", + width_cm = module_box_width_cm, + border_gp = grid::gpar(col = "black"), which = "row" ) right_anno <- ComplexHeatmap::HeatmapAnnotation( diff --git a/R/hc_utils.R b/R/hc_utils.R index dbfff51..2064233 100644 --- a/R/hc_utils.R +++ b/R/hc_utils.R @@ -2,14 +2,20 @@ #' @noRd get_cluster_colours <- function() { + # NOTE: colours must stay unique. Modules are coloured by index here, but + # `merge_clusters()` rebuilds its module table keyed by colour, so a repeated + # colour would silently collapse two distinct groups into one. The second + # "slategray" (position 39) was replaced with "purple" to keep the palette + # unique without shifting the indices of earlier colours. `unique()` guards + # against accidental duplicates being reintroduced later. col_vec <- c( "coral", "gold", "steelblue", "lightgreen", "turquoise", "plum", "maroon", "seagreen", "wheat", "slategray", "lightblue", "orchid", "darkgreen", "darkorange", "darkgrey", "indianred", "pink", "sandybrown", "khaki", "darkblue", "cadetblue", "greenyellow", "cyan", "thistle", "darkmagenta", "red", "blue", "green", "yellow", "brown", "black", "darkgoldenrod", - "cornsilk", "firebrick", "deeppink", "dodgerblue", "lightpink", "midnightblue", "slategray", "aquamarine", "chocolate", + "cornsilk", "firebrick", "deeppink", "dodgerblue", "lightpink", "midnightblue", "purple", "aquamarine", "chocolate", "darkred", "navy", "olivedrab", "peachpuff", "tomato", "snow" ) - return(col_vec) + return(base::unique(col_vec)) } .hc_match_axis_indices_with_duplicates <- function(axis_ids, requested_order) { @@ -3487,24 +3493,26 @@ hub_node_detection <- function(cluster, top, save, tree_layout, TF_only, plot) { return(hub_out) } - # vertex size (5 for hub, 3 for non-hub): - vertex_size <- base::lapply(igraph::V(g)$name, function(node) { - if (node %in% hub_out$hub_nodes) { - 5 - } else { - 3 - } - }) %>% base::unlist() + # vertex size: non-hubs small & uniform, hubs scaled by their combined + # centrality rank so the strongest hub reads as the largest node. + rank_lookup <- stats::setNames(hub_out$rank_df$sum, hub_out$rank_df$node) + hub_sums <- rank_lookup[hub_out$hub_nodes] + hub_sizes <- if (base::length(hub_sums) > 0 && + base::diff(base::range(hub_sums, na.rm = TRUE)) > 0) { + scales::rescale(hub_sums, to = base::c(6, 14)) + } else { + base::rep(9, base::length(hub_sums)) + } + base::names(hub_sizes) <- hub_out$hub_nodes + vertex_size <- base::vapply(igraph::V(g)$name, function(node) { + if (node %in% hub_out$hub_nodes) hub_sizes[[node]] else 2.5 + }, FUN.VALUE = base::numeric(1)) # layout: if (cluster == "all") { l <- hcobject[["integrated_output"]][["cluster_calc"]][["layout"]] } else { - if (tree_layout) { - l <- igraph::layout_as_tree(g) - } else { - l <- igraph::layout.lgl(g) # issues with >1 graph components? - } + l <- .hc_hub_network_layout(g, tree_layout = tree_layout) base::rownames(l) <- igraph::V(g)$name } @@ -3513,17 +3521,43 @@ hub_node_detection <- function(cluster, top, save, tree_layout, TF_only, plot) { igraph::V(g)$size <- vertex_size igraph::V(g)$label <- NA igraph::V(g)$color <- hub_out$colour_df$colour - # plot network: - network_with_labels( - network = g, - gene_labels = hub_out$hub_nodes, - gene_ranks = base::seq_along(hub_out$hub_nodes), - l = l, - label_offset = 10, - title = base::c(cluster, top), - save = save, - plot = plot - ) + # fade the non-hub background so the (fully opaque) hub nodes stand out. + is_hub <- igraph::V(g)$name %in% hub_out$hub_nodes + if (base::any(!is_hub)) { + igraph::V(g)$color[!is_hub] <- grDevices::adjustcolor( + igraph::V(g)$color[!is_hub], + alpha.f = 0.5 + ) + } + # plot network: modern ggplot2 style by default, with a graceful fallback to + # the classic base-graphics renderer when ggrepel is unavailable or the user + # opts out via options(hcocena.hub_network_style = "classic"). + use_modern <- !base::identical( + base::getOption("hcocena.hub_network_style", "modern"), "classic" + ) && .hc_has_modern_hub_plot_deps() + if (use_modern) { + .hc_hub_network_modern_plot( + network = g, + hub_nodes = hub_out$hub_nodes, + gene_ranks = base::seq_along(hub_out$hub_nodes), + layout = l, + centrality = stats::setNames(hub_out$rank_df$sum, hub_out$rank_df$node), + title = base::c(cluster, top), + save = save, + plot = plot + ) + } else { + network_with_labels( + network = g, + gene_labels = hub_out$hub_nodes, + gene_ranks = base::seq_along(hub_out$hub_nodes), + l = l, + label_offset = 10, + title = base::c(cluster, top), + save = save, + plot = plot + ) + } return(hub_out) @@ -3655,6 +3689,65 @@ combined_centrality <- function(network) { } +#' Hub Network Layout +#' +#' Layout for the per-cluster hub networks. Prefers graphlayouts' stress layout +#' (deterministic, cleanly separates components, avoids the hairball/cutoff +#' issues of `layout.lgl`), falling back to a weighted Fruchterman-Reingold and +#' finally `layout.lgl` so behaviour degrades gracefully when graphlayouts is +#' unavailable or a layout cannot be computed. +#' @noRd + +.hc_hub_network_layout <- function(g, tree_layout = FALSE) { + if (base::isTRUE(tree_layout)) { + return(igraph::layout_as_tree(g)) + } + weights <- .hc_graph_edge_weights(g) + l <- .hc_with_seed(1L, tryCatch( + { + if (base::requireNamespace("graphlayouts", quietly = TRUE)) { + graphlayouts::layout_with_stress(g) + } else { + igraph::layout_with_fr(g, weights = weights) + } + }, + error = function(e) NULL + )) + if (base::is.null(l) || !base::is.matrix(l) || base::nrow(l) != igraph::vcount(g)) { + l <- igraph::layout.lgl(g) + } + l +} + + +#' Rescale A Network Layout To A Stable Coordinate Range +#' +#' Maps layout coordinates into a `[0, target]` box while preserving the aspect +#' ratio (the larger of the two spans is scaled to `target`). This decouples +#' downstream label placement -- which uses an absolute offset -- from the +#' native coordinate scale of the layout backend, so stress/FR layouts are no +#' longer squished into a thin strip the way a fixed offset against their small +#' coordinate range would cause. +#' @noRd + +.hc_rescale_layout <- function(l, target = 100) { + l <- base::as.matrix(l) + if (base::nrow(l) == 0) { + return(l) + } + rng_x <- base::range(l[, 1], na.rm = TRUE) + rng_y <- base::range(l[, 2], na.rm = TRUE) + span <- base::max(base::diff(rng_x), base::diff(rng_y)) + if (!base::is.finite(span) || span <= 0) { + return(l) + } + scale <- target / span + l[, 1] <- (l[, 1] - rng_x[1]) * scale + l[, 2] <- (l[, 2] - rng_y[1]) * scale + l +} + + #' Weighted Degree Centrality #' #' Caclulates the weighted degree centrality of a node. Modified to be weighted from https://doi.org/10.1155/2019/9728742. @@ -3750,7 +3843,7 @@ weighted_BC <- function(network) { #' @noRd centrality_colours <- function(rank_df, network) { - mypalette <- grDevices::colorRampPalette(base::c("#FFF7EC", "#FC8D59", "#7F0000")) + mypalette <- grDevices::colorRampPalette(base::c("#FEE8C8", "#FC8D59", "#B30000")) colour_df <- base::data.frame( sum = rank_df$sum, id = rank_df$id, @@ -3761,6 +3854,122 @@ centrality_colours <- function(rank_df, network) { } +#' Modern Hub-Network Renderer Dependencies +#' +#' TRUE when the optional packages for the ggplot2 hub-network style are present. +#' @noRd + +.hc_has_modern_hub_plot_deps <- function() { + base::requireNamespace("ggplot2", quietly = TRUE) && + base::requireNamespace("ggrepel", quietly = TRUE) +} + + +#' Plot A Hub Network (Modern ggplot2 Style) +#' +#' Renders the per-cluster hub network as a clean ggplot2 figure: the full +#' network sits faintly in the background, the hub genes are drawn as bright +#' nodes sized and coloured by their combined centrality, and only the hubs are +#' labelled (directly at the node, with a white halo and short repelled +#' connectors via ggrepel). `coord_equal()` keeps the network undistorted. +#' Falls back to [network_with_labels()] is handled by the caller. +#' @noRd + +.hc_hub_network_modern_plot <- function(network, + hub_nodes, + gene_ranks, + layout, + centrality, + title, + save, + plot, + label_fontsize = 3.6, + node_size_range = base::c(3, 8)) { + plot_title <- base::paste0("cluster '", title[1], "' hub genes [top ", title[2], "]") + + node_names <- igraph::V(network)$name + coords <- base::as.matrix(layout) + if (!base::is.null(base::rownames(coords)) && + base::all(node_names %in% base::rownames(coords))) { + coords <- coords[node_names, , drop = FALSE] + } + + is_hub <- node_names %in% hub_nodes + nodes_df <- base::data.frame( + x = coords[, 1], + y = coords[, 2], + cent = base::as.numeric(centrality[node_names]), + is_hub = is_hub, + label = node_names, + stringsAsFactors = FALSE + ) + hubs_df <- nodes_df[nodes_df$is_hub, , drop = FALSE] + bg_df <- nodes_df[!nodes_df$is_hub, , drop = FALSE] + + el <- igraph::as_edgelist(network, names = FALSE) + edges_df <- if (base::nrow(el) > 0) { + base::data.frame( + x = coords[el[, 1], 1], y = coords[el[, 1], 2], + xend = coords[el[, 2], 1], yend = coords[el[, 2], 2] + ) + } else { + base::data.frame(x = base::numeric(0), y = base::numeric(0), + xend = base::numeric(0), yend = base::numeric(0)) + } + + p <- ggplot2::ggplot() + if (base::nrow(edges_df) > 0) { + p <- p + ggplot2::geom_segment( + data = edges_df, + ggplot2::aes(x = x, y = y, xend = xend, yend = yend), + colour = "grey80", alpha = 0.18, linewidth = 0.15 + ) + } + if (base::nrow(bg_df) > 0) { + p <- p + ggplot2::geom_point( + data = bg_df, ggplot2::aes(x = x, y = y), + colour = "grey75", alpha = 0.45, size = 0.9 + ) + } + p <- p + + ggplot2::geom_point( + data = hubs_df, + ggplot2::aes(x = x, y = y, fill = cent, size = cent), + shape = 21, colour = "white", stroke = 0.9 + ) + + ggplot2::scale_fill_viridis_c(option = "rocket", direction = -1, end = 0.92, name = "centrality") + + ggplot2::scale_size(range = node_size_range, guide = "none") + + ggrepel::geom_text_repel( + data = hubs_df, + ggplot2::aes(x = x, y = y, label = label), + fontface = "bold", size = label_fontsize, colour = "grey15", + bg.color = "white", bg.r = 0.18, + box.padding = 0.6, point.padding = 0.3, min.segment.length = 0, + segment.colour = "grey55", segment.size = 0.3, + max.overlaps = Inf, seed = 7 + ) + + ggplot2::coord_equal(clip = "off") + + ggplot2::labs(title = plot_title) + + ggplot2::theme_void(base_size = 13) + + ggplot2::theme( + plot.title = ggplot2::element_text(hjust = 0.5, face = "bold", size = 15), + plot.margin = ggplot2::margin(12, 18, 12, 18), + legend.position = "right" + ) + + if (isTRUE(save)) { + .hc_ggsave_pdf_png( + filename = .hc_output_file(base::paste0("Hub_genes_", title[1], "_module_network.pdf")), + plot = p, width = 10, height = 9, bg = "white" + ) + } + if (isTRUE(plot)) { + base::print(p) + } + base::invisible(p) +} + + #' Plots A Network With Labels #' #' Subroutine to find_hubs(). @@ -3769,12 +3978,23 @@ centrality_colours <- function(rank_df, network) { network_with_labels <- function(network, gene_labels, gene_ranks, l, label_offset, title, save, plot) { plot_title <- base::paste0("cluster '", title[1], "' hub genes [top ", title[2], "]") + # Normalise the layout to a stable coordinate scale and derive the side-label + # offset from the node spread. A fixed absolute offset squished small-scale + # layouts (stress/FR) into a thin vertical strip because it dwarfed the actual + # network extent; a proportional offset keeps a sensible aspect ratio for any + # layout backend. + l <- .hc_rescale_layout(l) + label_offset <- 0.55 * base::diff(base::range(l[, 1])) + new_genes_df <- base::data.frame(name = base::paste0("label_", base::seq_along(gene_labels)), label = gene_labels %>% base::as.character()) - mean_coord_x <- stats::median(l[, 1]) - mean_coord_y <- stats::median(l[, 2]) new_indeces <- base::match(new_genes_df$label, igraph::get.vertex.attribute(network)$name) new_genes_df$coords_x <- l[new_indeces, 1] new_genes_df$coords_y <- l[new_indeces, 2] + # Split the labels into the left/right (and up/down) columns around the median + # of the *hub* positions instead of all nodes, so the two columns stay roughly + # balanced even when the hubs cluster on one side of the layout. + mean_coord_x <- stats::median(new_genes_df$coords_x, na.rm = TRUE) + mean_coord_y <- stats::median(new_genes_df$coords_y, na.rm = TRUE) left_up <- dplyr::filter(new_genes_df, coords_x <= mean_coord_x & coords_y >= mean_coord_y) %>% dplyr::pull(., label) left_down <- dplyr::filter(new_genes_df, coords_x <= mean_coord_x & coords_y < mean_coord_y) %>% dplyr::pull(., label) right_up <- dplyr::filter(new_genes_df, coords_x > mean_coord_x & coords_y >= mean_coord_y) %>% dplyr::pull(., label) @@ -3847,6 +4067,21 @@ network_with_labels <- function(network, gene_labels, gene_ranks, l, label_offse } }) + # edge aesthetics: thin straight leader lines to the labels, co-expression + # edges scaled by weight and gently curved to reduce overplotting. + el <- igraph::get.edgelist(network2) + is_leader_edge <- el[, 2] %in% base::as.character(new_genes_df$name) + edge_weights <- .hc_graph_edge_weights(network2, default = NA_real_) + net_w <- edge_weights[!is_leader_edge] + new_edge_width <- base::rep(0.6, base::nrow(el)) + if (base::any(!is_leader_edge) && + base::diff(base::range(net_w, na.rm = TRUE)) > 0) { + new_edge_width[!is_leader_edge] <- scales::rescale(net_w, to = base::c(0.3, 2.2)) + } else if (base::any(!is_leader_edge)) { + new_edge_width[!is_leader_edge] <- 0.8 + } + new_edge_curved <- base::ifelse(is_leader_edge, 0, 0.12) + vertex_shape <- base::lapply(igraph::V(network2)$name, function(x) { if (x %in% igraph::V(network)$name) { "circle" @@ -3879,6 +4114,42 @@ network_with_labels <- function(network, gene_labels, gene_ranks, l, label_offse 0 } }) %>% base::unlist() + + # outline the hub nodes (the labelled genes) so they stand out from the faded + # background; leave the rest borderless to avoid speckling dense networks. + is_hub_node <- igraph::V(network2)$name %in% base::as.character(gene_labels) + vertex_frame_color <- base::ifelse(is_hub_node, "black", NA) + + # Normalise l2 into a centred square and plot with rescale = FALSE / asp = 1. + # igraph's default per-axis rescale stretches each axis to [-1, 1] + # independently; because the side labels widen only the x-range, that made + # the network render roughly twice as tall as wide. Scaling both axes by the + # same factor preserves the network's true aspect ratio. + # Centre the frame on the network-node centroid (not the bounding-box midpoint) + # so the dense core sits in the middle and a few peripheral nodes no longer pull + # the composition off to one side. Scale both axes by the same factor to keep + # the network undistorted (plotted with rescale = FALSE / asp = 1 below). + net_rows <- base::seq_len(igraph::vcount(network)) + l2_cx <- stats::median(l2[net_rows, 1], na.rm = TRUE) + l2_cy <- stats::median(l2[net_rows, 2], na.rm = TRUE) + l2_scale <- base::max(base::diff(base::range(l2[, 1])), base::diff(base::range(l2[, 2]))) / 2 + if (base::is.finite(l2_scale) && l2_scale > 0) { + l2[, 1] <- (l2[, 1] - l2_cx) / l2_scale + l2[, 2] <- (l2[, 2] - l2_cy) / l2_scale + } + # window reaches the farthest label on each axis (extra horizontal room for the + # label text) so nothing is clipped while the core stays centred. + x_reach <- base::max(base::abs(l2[, 1]), na.rm = TRUE) + y_reach <- base::max(base::abs(l2[, 2]), na.rm = TRUE) + plot_xlim <- base::c(-x_reach, x_reach) * 1.18 + plot_ylim <- base::c(-y_reach, y_reach) * 1.08 + if (!base::all(base::is.finite(plot_xlim)) || base::diff(plot_xlim) <= 0) { + plot_xlim <- base::c(-1.2, 1.2) + } + if (!base::all(base::is.finite(plot_ylim)) || base::diff(plot_ylim) <= 0) { + plot_ylim <- base::c(-1.2, 1.2) + } + if (save == TRUE) { .hc_export_single_page_plot( file = .hc_output_file(base::paste0("Hub_genes_", title[1], "_module_network.pdf")), @@ -3891,9 +4162,16 @@ network_with_labels <- function(network, gene_labels, gene_ranks, l, label_offse vertex.label = new_labels, vertex.label.cex = 1.5, layout = l2, + rescale = FALSE, + xlim = plot_xlim, + ylim = plot_ylim, + asp = 1, vertex.label.dist = 1, vertex.shape = vertex_shape, + vertex.frame.color = vertex_frame_color, edge.color = new_edge_color, + edge.width = new_edge_width, + edge.curved = new_edge_curved, vertex.label.color = new_label_color ) graphics::title(plot_title, cex.main = 3) @@ -3904,7 +4182,10 @@ network_with_labels <- function(network, gene_labels, gene_ranks, l, label_offse igraph::plot.igraph(network2, vertex.size = vertex_size, vertex.label = new_labels, vertex.label.cex = 0.75, layout = l2, main = plot_title, vertex.label.dist = 1, vertex.shape = vertex_shape, - edge.color = new_edge_color, vertex.label.color = new_label_color + rescale = FALSE, xlim = plot_xlim, ylim = plot_ylim, asp = 1, + vertex.frame.color = vertex_frame_color, edge.color = new_edge_color, + edge.width = new_edge_width, edge.curved = new_edge_curved, + vertex.label.color = new_label_color ) } } diff --git a/README.md b/README.md index 93e721a..88fad53 100644 --- a/README.md +++ b/README.md @@ -76,6 +76,8 @@ The container prepares a workspace at `/home/rstudio/hcocena` and includes: - visible workflow notebooks under `/home/rstudio/hcocena/github_workflows/` including `hcocena_main.Rmd` and `hcocena_satellite.Rmd` - preinstalled optional packages for common workflows, including + `DESeq2`, `limma`, `sva`, `edgeR`, `tximport`, `apeglm`, `ashr`, + `EnhancedVolcano`, `pheatmap`, `org.Hs.eg.db`, `org.Mm.eg.db`, `CALIBERrfimpute`, `RCy3`, `SpatialExperiment`, and `GSVA` - empty `count_data`, `annotation_data`, and `output` directories @@ -120,6 +122,44 @@ browseVignettes("hcocena") The package ships toy data and prepared example objects in `inst/extdata` to support documentation, testing, and manual smoke tests. +## Real-data regression checks + +Local real-data checks live behind an explicit opt-in runner. By default, the +runner looks for the STAR protocol data in `../data` relative to this repository +and writes ignored outputs to `realdata-output/`. + +`quick` uses the real Array/RNA-seq data with the top 2000 variable genes per +layer. It runs integration, clustering, module heatmaps, module splitting, +Hallmark/KEGG enrichment, and a module-label font-size probe. `full` keeps the +same checks but uses the full imported gene set and broader enrichment defaults. + +```bash +Rscript scripts/run_realdata_regression.R --mode quick +Rscript scripts/run_realdata_regression.R --mode quick --update-reference +Rscript scripts/run_realdata_regression.R --mode full +``` + +Set `HCOCENA_REALDATA_DIR` to point at another data directory. Normal package +tests skip the real-data run; enable it explicitly with: + +```bash +HCOCENA_RUN_REALDATA=true HCOCENA_REALDATA_MODE=quick Rscript -e "testthat::test_file('tests/testthat/test-realdata-regression.R')" +``` + +PowerShell: + +```powershell +$env:HCOCENA_RUN_REALDATA = "true"; $env:HCOCENA_REALDATA_MODE = "quick"; Rscript -e "testthat::test_file('tests/testthat/test-realdata-regression.R')" +``` + +The exported CSV/JSON artifacts are designed as golden-master inputs for later +R/Python comparisons: QC metrics, edge lists, GFC matrices, module tables, +split diagnostics, enrichment tables, heatmap matrices, and a module-label +font-size probe. Each run also writes a `visual_check_report_.pdf` with +the generated plot variants, non-default plot parameters, and image pages for +manual inspection, including the pre-split heatmap, post-split heatmap, and +combined enrichment plot. + ## GitHub workflows This repository also ships longer GitHub-oriented walkthroughs: diff --git a/docker/DOCKERHUB_OVERVIEW.md b/docker/DOCKERHUB_OVERVIEW.md index 5d53680..7a1184d 100644 --- a/docker/DOCKERHUB_OVERVIEW.md +++ b/docker/DOCKERHUB_OVERVIEW.md @@ -13,14 +13,19 @@ integration and downstream analysis of transcriptomics datasets. - `02_hcocena_satellite.Rmd` - `03_hcocena_main_seq_only.Rmd` - additional notebooks under `/home/rstudio/hcocena/github_workflows/` +- RNA-seq and differential-expression packages including `DESeq2`, `limma`, + `sva`, `edgeR`, `tximport`, `apeglm`, `ashr`, `EnhancedVolcano`, + `pheatmap`, `org.Hs.eg.db`, and `org.Mm.eg.db` ## Recommended tags - `latest` for the current image and quick-start commands -- `1.97` for a pinned, reproducible setup +- `1.99` for a pinned, reproducible setup Older tags still kept for older reproducible runs: +- `1.98` +- `1.97` - `1.96` - `1.95` - `1.94` @@ -42,7 +47,7 @@ Run RStudio Server: docker run --rm -p 8787:8787 -e PASSWORD=hcocena therealtomek/hcocena:latest ``` -For reproducible runs, replace `latest` with a pinned tag such as `1.97`. +For reproducible runs, replace `latest` with a pinned tag such as `1.99`. Then open: diff --git a/docker/Dockerfile b/docker/Dockerfile index fc5760e..e64b3dc 100644 --- a/docker/Dockerfile +++ b/docker/Dockerfile @@ -36,12 +36,14 @@ RUN apt-get update && apt-get install -y --no-install-recommends \ && rm -rf /var/lib/apt/lists/* WORKDIR /opt/hcocena -COPY . /opt/hcocena RUN R -q -e "options(timeout = 1200); install.packages(c('remotes', 'BiocManager', 'graphlayouts', 'ellmer'), repos = 'https://cloud.r-project.org')" && \ - R -q -e "options(timeout = 1200); repos <- BiocManager::repositories(); repos['CRAN'] <- 'https://cloud.r-project.org'; options(repos = repos); pkgs <- c('ape', 'treeio', 'ggtree', 'GO.db', 'GOSemSim', 'DOSE', 'enrichplot', 'clusterProfiler', 'CALIBERrfimpute', 'decoupleR', 'dorothea', 'mice', 'progeny', 'RCy3', 'SpatialExperiment', 'GSVA'); for (attempt in seq_len(3)) { missing <- pkgs[!vapply(pkgs, requireNamespace, logical(1), quietly = TRUE)]; if (!length(missing)) break; message('Docker dependency install attempt ', attempt, ': ', paste(missing, collapse = ', ')); BiocManager::install(missing, ask = FALSE, update = FALSE) }; status <- vapply(pkgs, requireNamespace, logical(1), quietly = TRUE); print(status); if (!all(status)) stop('Missing Docker dependency packages: ', paste(names(status)[!status], collapse = ', '))" && \ - R -q -e "options(timeout = 1200); repos <- BiocManager::repositories(); repos['CRAN'] <- 'https://cloud.r-project.org'; options(repos = repos); remotes::install_local('/opt/hcocena', dependencies = NA, upgrade = 'never')" && \ - R -q -e "status <- vapply(c('graphlayouts', 'ellmer', 'clusterProfiler', 'CALIBERrfimpute', 'decoupleR', 'dorothea', 'mice', 'progeny', 'RCy3', 'SpatialExperiment', 'GSVA', 'hcocena'), requireNamespace, logical(1), quietly = TRUE); print(status); if (!all(status)) stop('Missing Docker packages: ', paste(names(status)[!status], collapse = ', '))" && \ + R -q -e "options(timeout = 1200); repos <- BiocManager::repositories(); repos['CRAN'] <- 'https://cloud.r-project.org'; options(repos = repos); pkgs <- c('ape', 'treeio', 'ggtree', 'GO.db', 'org.Hs.eg.db', 'org.Mm.eg.db', 'GOSemSim', 'DOSE', 'enrichplot', 'clusterProfiler', 'DESeq2', 'limma', 'sva', 'edgeR', 'tximport', 'apeglm', 'ashr', 'EnhancedVolcano', 'pheatmap', 'CALIBERrfimpute', 'decoupleR', 'dorothea', 'mice', 'progeny', 'RCy3', 'SpatialExperiment', 'GSVA'); for (attempt in seq_len(3)) { missing <- pkgs[!vapply(pkgs, requireNamespace, logical(1), quietly = TRUE)]; if (!length(missing)) break; message('Docker dependency install attempt ', attempt, ': ', paste(missing, collapse = ', ')); BiocManager::install(missing, ask = FALSE, update = FALSE) }; status <- vapply(pkgs, requireNamespace, logical(1), quietly = TRUE); print(status); if (!all(status)) stop('Missing Docker dependency packages: ', paste(names(status)[!status], collapse = ', '))" + +COPY . /opt/hcocena + +RUN R -q -e "options(timeout = 1200); repos <- BiocManager::repositories(); repos['CRAN'] <- 'https://cloud.r-project.org'; options(repos = repos); remotes::install_local('/opt/hcocena', dependencies = NA, upgrade = 'never')" && \ + R -q -e "status <- vapply(c('graphlayouts', 'ellmer', 'clusterProfiler', 'DESeq2', 'limma', 'sva', 'edgeR', 'tximport', 'apeglm', 'ashr', 'EnhancedVolcano', 'pheatmap', 'org.Hs.eg.db', 'org.Mm.eg.db', 'CALIBERrfimpute', 'decoupleR', 'dorothea', 'mice', 'progeny', 'RCy3', 'SpatialExperiment', 'GSVA', 'hcocena'), requireNamespace, logical(1), quietly = TRUE); print(status); if (!all(status)) stop('Missing Docker packages: ', paste(names(status)[!status], collapse = ', '))" && \ R -q -e "packageVersion('hcocena')" RUN mkdir -p /home/rstudio/count_data && \ diff --git a/docker/README.md b/docker/README.md index cf9b275..81b5b7a 100644 --- a/docker/README.md +++ b/docker/README.md @@ -18,7 +18,9 @@ The image includes: - workflow notebooks directly in the workspace root: `01_hcocena_main.Rmd`, `02_hcocena_satellite.Rmd`, `03_hcocena_main_seq_only.Rmd` - preinstalled optional packages for common workflows, including - `graphlayouts`, `ellmer`, `CALIBERrfimpute`, `RCy3`, `SpatialExperiment`, and `GSVA` + `graphlayouts`, `ellmer`, `DESeq2`, `limma`, `sva`, `edgeR`, `tximport`, + `apeglm`, `ashr`, `EnhancedVolcano`, `pheatmap`, `org.Hs.eg.db`, + `org.Mm.eg.db`, `CALIBERrfimpute`, `RCy3`, `SpatialExperiment`, and `GSVA` - empty `count_data`, `annotation_data`, and `output` directories RStudio Server is exposed on port `8787` by the base image. diff --git a/docker/project_workflows/hcocena_satellite.Rmd b/docker/project_workflows/hcocena_satellite.Rmd index dc58d9b..5c8f61e 100644 --- a/docker/project_workflows/hcocena_satellite.Rmd +++ b/docker/project_workflows/hcocena_satellite.Rmd @@ -381,88 +381,88 @@ hc <- hc_plot_cluster_heatmap( -Purpose: summarize module biology using Gemini, OpenAI, or vLLM. +Purpose: summarize module biology with an LLM. One generic function, +`hc_module_function_llm()`, drives every provider — pick the model via +`provider`/`model`, give `base_url` for a local server, and `api_key` when needed. Key parameters: - `module`: one label (e.g. `"M1"`) or `"all"`. - `biological_context`: study context to guide interpretation. - `system_instruction`: prompt constraints. -- Provider/model params: - `llm`, `gemini_model`, `openai_model`, `vllm_model`, `vllm_base_url`. +- `provider`: `"gemini"`, `"claude"`, `"openai"`, or `"vllm"` (local server). +- `model` / `base_url` / `api_key`: model code, local endpoint, API key. + +List the models each provider currently offers: + +```{r eval=FALSE, message=TRUE, warning=FALSE} +hc_list_llm_models(c("gemini", "openai", "claude")) +hc_list_llm_models("vllm", base_url = "http://bn-hl-compute-8mi210-03.dzne.de:8000/v1") +``` ```{r eval=FALSE, fig.width=13, fig.height=10, out.width="100%", message=TRUE, warning=FALSE} biological_context <- "maturation of monocytes in preterm infants over the first year of life" -gemini_api_key <- Sys.getenv("GEMINI_API_KEY") -openai_api_key <- Sys.getenv("OPENAI_API_KEY") -anthropic_api_key <- Sys.getenv("ANTHROPIC_API_KEY") -system_instruction <- paste("You are an expert transcriptomics analyst.", "Task: Infer the overarching biological program from the provided list of co-expressed genes.", "Use the provided biological context to inform your interpretation and synthesize it into a precise biological description.", "For the 'contextual_state' field, strictly describe the specific cellular or transcriptional state.", "Constraint 1: Start your descriptions directly with the relevant nouns or biological terms. Completely omit filler phrases like 'this module' or 'the genes'.", "Constraint 2: Formulate your answer as a direct biological interpretation based on the supplied genes, without claiming any statistical enrichment tests.", "Output: Return ONLY a valid JSON object matching the requested schema. Do not output markdown formatting, code blocks, or any conversational text.") +gemini_api_key <- Sys.getenv("GEMINI_API_KEY") +openai_api_key <- Sys.getenv("OPENAI_API_KEY") +anthropic_api_key <- Sys.getenv("ANTHROPIC_API_KEY") -hc <- hc_module_function_llm( - hc, - module = "M1", - biological_context = biological_context, - system_instruction = system_instruction, - llm = "gemini", - gemini_model = "gemini-3.1-pro-preview", - api_key = gemini_api_key, - timeout_sec = NULL, - pause_sec = NULL, - continue_on_error = TRUE +system_instruction <- paste( + "You are an expert transcriptomics analyst.", + "Task: Infer the overarching biological program from the provided list of co-expressed genes.", + "Use the provided biological context to inform your interpretation and synthesize it into a precise biological description.", + "For the 'contextual_state' field, strictly describe the specific cellular or transcriptional state.", + "Constraint 1: Start your descriptions directly with the relevant nouns or biological terms. Completely omit filler phrases like 'this module' or 'the genes'.", + "Constraint 2: Formulate your answer as a direct biological interpretation based on the supplied genes, without claiming any statistical enrichment tests.", + "Output: Return ONLY a valid JSON object matching the requested schema. Do not output markdown formatting, code blocks, or any conversational text." ) -hc <- hc_module_function_llm( - hc, - module = "all", - biological_context = biological_context, - system_instruction = system_instruction, - llm = "gemini", - gemini_model = "gemini-3.1-pro-preview", - api_key = gemini_api_key, - timeout_sec = NULL, - pause_sec = NULL, - continue_on_error = TRUE -) -# -# #Optional providers: -# hc <- hc_module_function_llm(hc, module = "all", llm = "openai", -# openai_model = "gpt-4.1-mini", api_key = openai_api_key) +## Configure one entry per model: provider + model, plus base_url (for a local +## server) and api_key where needed. Comment lines in/out to (de)activate models. +models <- list( + gemini = list(provider = "gemini", model = "gemini-3.1-pro-preview", + api_key = gemini_api_key), -hc <- hc_module_function_vllm( - hc, - module = "all", - biological_context = biological_context, - system_instruction = system_instruction, - vllm_model = "/data/qwen3.5-35b-a3b", - vllm_base_url = "http://bn-hl-compute-8mi210-03.dzne.de:8000/v1", - api_key = "EMPTY", - temperature = 0.8, - timeout_sec = 600, - max_genes = 1250, - continue_on_error = TRUE -) + qwen = list(provider = "vllm", model = "/data/qwen3.5-35b-a3b", + base_url = "http://bn-hl-compute-8mi210-03.dzne.de:8000/v1", + api_key = "EMPTY", + temperature = 0.8, timeout_sec = 600, max_genes = 1250), -hc <- hc_module_function_llm( - hc = hc, - module = "all", - biological_context = biological_context, - system_instruction = system_instruction, - llm = "claude", - claude_model = "claude-sonnet-4-6", - api_key = anthropic_api_key, - timeout_sec = NULL, - pause_sec = NULL, - continue_on_error = TRUE -) + claude = list(provider = "claude", model = "claude-opus-4-6", + api_key = anthropic_api_key) + # openai = list(provider = "openai", model = "gpt-4.1-mini", + # api_key = openai_api_key) +) +## One generic call per model; each result is stored in its own satellite slot. +for (nm in names(models)) { + message("== LLM module function: ", nm, " ==") + hc <- tryCatch( + do.call(hc_module_function_llm, c( + list( + hc = hc, + module = "all", + biological_context = biological_context, + system_instruction = system_instruction, + slot_name = paste0("llm_module_function_", nm), + continue_on_error = TRUE + ), + models[[nm]] + )), + error = function(e) { + warning("Skipping '", nm, "': ", conditionMessage(e), call. = FALSE) + hc + } + ) +} -plots <- hc_plot_module_function_llm(hc) +## Plot one model's interpretation (slot = llm_module_function_). +plots <- hc_plot_module_function_llm(hc, slot_name = "llm_module_function_gemini") print(plots$general_processes) print(plots$contextual_state) print(plots$key_regulators) -#View(hc@satellite$llm_module_function_summary) +# View(hc@satellite$llm_module_function_gemini_summary) ``` diff --git a/man/hc_module_function_llm.Rd b/man/hc_module_function_llm.Rd index 8235301..a13a145 100644 --- a/man/hc_module_function_llm.Rd +++ b/man/hc_module_function_llm.Rd @@ -13,17 +13,31 @@ hc_module_function_llm( context = NULL, biological_context = NULL, label = NULL, + provider = NULL, llm = c("gemini", "claude", "openai", "chatgpt", "vllm"), api_key = NULL, gemini_model = NULL, claude_model = NULL, vllm_model = NULL, vllm_base_url = NULL, + base_url = NULL, model = NULL, max_genes = Inf, temperature = 0.2, timeout_sec = 60, pause_sec = 0, + use_rag = FALSE, + rag_query = NULL, + rag_url = NULL, + rag_category = "DoRAG", + rag_top_k = 10, + rag_themes = NULL, + rag_timeout_sec = 120, + rag_connect_timeout_sec = 30, + rag_min_relevance = NULL, + rag_max_context_chars = 12000, + rag_continue_on_error = FALSE, + compare_interpretation_levels = FALSE, continue_on_error = FALSE, save_to_hc = !base::is.null(hc), slot_name = "llm_module_function", @@ -56,6 +70,9 @@ context for the interpretation. Preferred over \code{context}.} module name or \code{"custom_geneset"}. Can only be used for a single module or a free gene set.} +\item{provider}{Provider-neutral alias for \code{llm}. When supplied it takes +precedence over \code{llm}. Accepts the same values as \code{llm}.} + \item{llm}{LLM provider. Supported values are \code{"gemini"}, \code{"claude"}, \code{"anthropic"} (alias for \code{"claude"}), \code{"openai"}, \code{"chatgpt"} (alias for \code{"openai"}), or \code{"vllm"} for a local OpenAI-compatible vLLM server.} @@ -77,6 +94,10 @@ Defaults to \code{"Qwen/Qwen2.5-VL-32B-Instruct"}.} \item{vllm_base_url}{Base URL for the local OpenAI-compatible vLLM server. Only used when \code{llm = "vllm"}. Defaults to \code{"http://localhost:8000/v1"}.} +\item{base_url}{Provider-neutral alias for the local server endpoint. When +supplied with \code{llm = "vllm"} it overrides \code{vllm_base_url}; for cloud +providers it is ignored (with a warning).} + \item{model}{Optional generic model code. Mainly useful for OpenAI, or as a provider-agnostic override.} @@ -92,6 +113,52 @@ requests unless overridden explicitly.} \item{pause_sec}{Pause in seconds between module requests. Useful for \code{module = "all"}. Use \code{0} or \code{NULL} for no pause. Default is \code{0}.} +\item{use_rag}{Logical. If \code{TRUE}, retrieves literature passages from the +DoRAG raw retrieval API and injects them into the LLM prompt as optional +supporting context. Default is \code{FALSE}.} + +\item{rag_query}{Optional retrieval query. If \code{NULL}, a query is built from +the biological context, module label, and submitted genes. A single string +is reused for all modules; a named character vector can map module labels +to queries; a character vector with one entry per module is used in order. +Advanced users may pass a function with arguments \code{label}, \code{module}, +\code{genes}, and \code{context_text}.} + +\item{rag_url}{DoRAG raw retrieval API endpoint. If \code{NULL}, uses +\code{HCOCENA_RAG_URL} or +\code{"https://limesbcnr-007901.iaas.uni-bonn.de/api/rag/query"}.} + +\item{rag_category}{Retrieval category sent to DoRAG. Default is \code{"DoRAG"}.} + +\item{rag_top_k}{Number of passages requested after reranking. Default is +\code{10}.} + +\item{rag_themes}{Optional biological theme filter passed to DoRAG, for +example \code{"Microbiome Science"} or \code{c("Microbiome Science", "Innate Immunity")}.} + +\item{rag_timeout_sec}{DoRAG request timeout in seconds. Use \code{0} or \code{NULL} +to disable the timeout. Default is \code{120}.} + +\item{rag_connect_timeout_sec}{DoRAG connection timeout in seconds. This +controls how long to wait while establishing the TCP connection to the +RAG server. Use \code{0} or \code{NULL} to use the system default. Default is \code{30}.} + +\item{rag_min_relevance}{Optional minimum rerank score. Retrieved passages +below this score are omitted from the prompt after retrieval.} + +\item{rag_max_context_chars}{Maximum number of characters of formatted DoRAG +context injected into each LLM prompt. Use \code{Inf} or \code{NULL} to disable this +limit. Default is \code{12000}.} + +\item{rag_continue_on_error}{Logical. If \code{TRUE}, a failed DoRAG request is +stored in the result and the LLM interpretation continues without RAG +passages. Default is \code{FALSE}.} + +\item{compare_interpretation_levels}{Logical. If \code{TRUE}, stores separate +LLM interpretations for comparison: \code{without_context}, \code{with_context}, +and, when \code{use_rag = TRUE}, \code{with_rag}. The top-level \code{response} remains +the final requested interpretation. Default is \code{FALSE}.} + \item{continue_on_error}{Logical. If \code{TRUE}, continue with the next module when one request fails and store the error in the result summary.} @@ -132,7 +199,6 @@ The function name is retained for backward compatibility, but it can now use Gemini, OpenAI-compatible APIs, and local vLLM servers via \code{llm}. } \examples{ -hc <- hc_init() \donttest{ if (nzchar(Sys.getenv("OPENAI_API_KEY"))) { res <- hc_module_function_llm( diff --git a/man/hc_plot_cluster_heatmap_view.Rd b/man/hc_plot_cluster_heatmap_view.Rd new file mode 100644 index 0000000..e017c13 --- /dev/null +++ b/man/hc_plot_cluster_heatmap_view.Rd @@ -0,0 +1,75 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/hc_plot_cluster_heatmap_view.R +\name{hc_plot_cluster_heatmap_view} +\alias{hc_plot_cluster_heatmap_view} +\title{Plot a fixed-scale view of selected modules from the module heatmap} +\usage{ +hc_plot_cluster_heatmap_view( + hc, + modules, + file_name = "Heatmap_modules_view.pdf", + col_order = NULL, + cluster_columns = FALSE, + cluster_rows = FALSE, + gfc_colors = NULL, + gfc_scale_limits = NULL, + preserve_full_heatmap = TRUE, + store_view = TRUE, + view_name = NULL, + ... +) +} +\arguments{ +\item{hc}{A \code{HCoCenaExperiment}.} + +\item{modules}{Character/numeric vector of modules to show. Accepts current +module labels such as \code{"M2"} or \code{"M2.1"}, module colors, or numeric module +indices in the current heatmap order.} + +\item{file_name}{Optional export file name for the view. Use \code{FALSE} to skip +file export.} + +\item{col_order}{Optional condition order. If \code{NULL}, the stored full +heatmap column order is reused when available.} + +\item{cluster_columns}{Logical. Whether to cluster columns in the view.} + +\item{cluster_rows}{Logical. Whether to cluster rows in the view. Defaults to +\code{FALSE} so the \code{modules} vector controls the displayed order.} + +\item{gfc_colors}{Optional GFC palette override. If \code{NULL}, the stored full +heatmap palette is reused.} + +\item{gfc_scale_limits}{Optional GFC scale override. If \code{NULL}, the stored +full heatmap scale is reused; if unavailable, a full-cache/global fallback +is used.} + +\item{preserve_full_heatmap}{Logical. If \code{TRUE}, restore the pre-existing +full heatmap cache in \code{hc} after writing the view.} + +\item{store_view}{Logical. If \code{TRUE}, store lightweight metadata for the +view under \code{hc@integration@cluster$heatmap_views}.} + +\item{view_name}{Optional name used when \code{store_view = TRUE}.} + +\item{...}{Additional plotting arguments forwarded to +\code{\link[=hc_plot_cluster_heatmap]{hc_plot_cluster_heatmap()}}.} +} +\value{ +Updated \code{HCoCenaExperiment}. +} +\description{ +Creates a module heatmap view for selected modules in the exact order +supplied by \code{modules}. The module labels and module-color annotations are +preserved, and the GFC expression color scale is inherited from the full +heatmap cache instead of being recomputed from the subset. +} +\examples{ +\donttest{ +hc <- hc_plot_cluster_heatmap_view( + hc, + modules = c("M2", "M3", "M4", "M8"), + file_name = "Heatmap_modules_selected.pdf" +) +} +} diff --git a/man/hc_plot_module_function_llm.Rd b/man/hc_plot_module_function_llm.Rd index f66baef..00122a8 100644 --- a/man/hc_plot_module_function_llm.Rd +++ b/man/hc_plot_module_function_llm.Rd @@ -42,8 +42,14 @@ hc_plot_module_function_gemini(...) \item{modules}{Optional character vector to subset modules.} \item{fields}{Character vector selecting which LLM fields to plot. Supported -values are \code{"general_processes"}, \code{"contextual_state"}, and -\code{"key_regulators"}. Defaults to all three.} +values include \code{"general_processes"}, \code{"contextual_state"}, +\code{"key_regulators"}, and comparison fields from +\code{compare_interpretation_levels = TRUE}, such as +\code{"without_context_contextual_state"}, \code{"with_context_contextual_state"}, +and \code{"with_rag_contextual_state"}. The aliases +\code{"contextual_state_without_context"}, \code{"contextual_state_with_context"}, +and \code{"contextual_state_rag"} are also accepted. Defaults to the three +top-level fields.} \item{max_chars}{Maximum number of characters shown per term. Default is \code{90}.} diff --git a/man/hc_split_modules.Rd b/man/hc_split_modules.Rd index 4ab4802..39a8307 100644 --- a/man/hc_split_modules.Rd +++ b/man/hc_split_modules.Rd @@ -33,7 +33,10 @@ One of \code{"cluster_leiden"} (default), \code{"cluster_louvain"}, \item{no_of_iterations}{Number of Leiden iterations (used only for Leiden).} -\item{resolution}{Leiden resolution (used only for Leiden).} +\item{resolution}{Leiden resolution (used only for Leiden). Use either one +positive value for all selected modules or one positive value per selected +module in the same order as \code{modules}. Named vectors may use module labels +or module colors.} \item{resolution_grid}{Optional numeric vector of candidate resolutions to test before splitting. For each candidate, hCoCena reports how many diff --git a/scripts/run_realdata_regression.R b/scripts/run_realdata_regression.R new file mode 100644 index 0000000..c81774f --- /dev/null +++ b/scripts/run_realdata_regression.R @@ -0,0 +1,1562 @@ +#!/usr/bin/env Rscript + +# Real-data regression runner for local hCoCena development. +# +# The runner is intentionally opt-in: real data and generated outputs are kept +# outside package tests unless HCOCENA_RUN_REALDATA is enabled. + +.hcr_bool <- function(x, default = FALSE) { + if (is.null(x) || length(x) == 0 || is.na(x) || !nzchar(x)) { + return(default) + } + tolower(trimws(as.character(x[[1]]))) %in% c("1", "true", "yes", "y", "on") +} + +.hcr_scalar <- function(x, default = NULL) { + if (is.null(x) || length(x) == 0 || is.na(x) || !nzchar(as.character(x[[1]]))) { + return(default) + } + as.character(x[[1]]) +} + +.hcr_parse_args <- function(args = commandArgs(trailingOnly = TRUE)) { + out <- list() + i <- 1L + while (i <= length(args)) { + arg <- args[[i]] + if (!startsWith(arg, "--")) { + stop("Unexpected positional argument: ", arg, call. = FALSE) + } + key <- sub("^--", "", arg) + value <- "true" + if (grepl("=", key, fixed = TRUE)) { + parts <- strsplit(key, "=", fixed = TRUE)[[1]] + key <- parts[[1]] + value <- paste(parts[-1], collapse = "=") + } else if (i < length(args) && !startsWith(args[[i + 1L]], "--")) { + i <- i + 1L + value <- args[[i]] + } + key <- gsub("-", "_", key, fixed = TRUE) + out[[key]] <- value + i <- i + 1L + } + out +} + +.hcr_repo_root <- function(start = getwd()) { + path <- normalizePath(start, winslash = "/", mustWork = TRUE) + repeat { + desc <- file.path(path, "DESCRIPTION") + if (file.exists(desc)) { + dcf <- tryCatch(read.dcf(desc), error = function(e) NULL) + if (!is.null(dcf) && identical(unname(dcf[1, "Package"]), "hcocena")) { + return(path) + } + } + parent <- dirname(path) + if (identical(parent, path)) { + stop("Could not find the hcocena repository root from: ", start, call. = FALSE) + } + path <- parent + } +} + +.hcr_slash <- function(path) { + path <- normalizePath(path, winslash = "/", mustWork = FALSE) + if (grepl("/$", path)) path else paste0(path, "/") +} + +.hcr_default_data_dir <- function(repo_root) { + env <- Sys.getenv("HCOCENA_REALDATA_DIR", unset = "") + if (nzchar(env)) { + return(env) + } + file.path(dirname(repo_root), "data") +} + +.hcr_default_reference_dir <- function(repo_root, mode) { + env <- Sys.getenv("HCOCENA_REALDATA_REFERENCE_DIR", unset = "") + if (nzchar(env)) { + return(file.path(env, mode)) + } + file.path(repo_root, "realdata-reference", mode) +} + +.hcr_load_package <- function(repo_root, load_source = TRUE) { + if (isTRUE(load_source) && requireNamespace("pkgload", quietly = TRUE)) { + pkgload::load_all(repo_root, quiet = TRUE) + return(invisible("source")) + } + suppressPackageStartupMessages(library(hcocena)) + invisible("installed") +} + +.hcr_mode_config <- function(mode) { + mode <- match.arg(mode, c("quick", "full")) + if (identical(mode, "quick")) { + return(list( + mode = mode, + top_var = c(2000, 2000), + min_corr = c(0.90, 0.90), + range_cutoff_length = c(5, 5), + cutoffs = c(0.90, 0.90), + min_nodes_number_for_network = 5, + min_nodes_number_for_cluster = 5, + cluster_algo = "cluster_fast_greedy", + no_of_iterations = 1, + resolution = 0.1, + read_supplementary = TRUE, + run_enrichment = TRUE, + enrichment_gene_sets = c("Hallmark", "Kegg"), + enrichment_top = 3, + run_module_split = TRUE, + split_cluster_algo = "cluster_fast_greedy", + split_resolution = 0.1, + split_resolution_grid = c(0.05, 0.1, 0.2), + split_min_submodule_size = 5, + plot_layer_heatmaps = FALSE, + module_label_fontsize = 11, + module_label_pt_size = 0.55, + module_box_width_cm = 0.80 + )) + } + + list( + mode = mode, + top_var = c("all", "all"), + min_corr = c(0.90, 0.90), + range_cutoff_length = c(100, 100), + cutoffs = c(0.90, 0.90), + min_nodes_number_for_network = 15, + min_nodes_number_for_cluster = 15, + cluster_algo = "cluster_leiden", + no_of_iterations = 2, + resolution = 0.1, + read_supplementary = TRUE, + run_enrichment = .hcr_bool(Sys.getenv("HCOCENA_REALDATA_FULL_ENRICHMENT", unset = "true"), TRUE), + enrichment_gene_sets = c("Hallmark", "Kegg", "Go"), + enrichment_top = 5, + run_module_split = TRUE, + split_cluster_algo = "cluster_fast_greedy", + split_resolution = 0.1, + split_resolution_grid = c(0.05, 0.1, 0.2), + split_min_submodule_size = 15, + plot_layer_heatmaps = FALSE, + module_label_fontsize = 11, + module_label_pt_size = 0.55, + module_box_width_cm = 0.80 + ) +} + +.hcr_check_inputs <- function(data_dir, reference_dir) { + data_dir <- normalizePath(data_dir, winslash = "/", mustWork = TRUE) + required <- c( + "data_seq_processed.txt", + "annotation_seq.txt", + "data_array_processed.txt", + "annotation_array.txt" + ) + missing <- required[!file.exists(file.path(data_dir, required))] + if (length(missing) > 0) { + stop( + "Real-data directory is missing required files: ", + paste(missing, collapse = ", "), + call. = FALSE + ) + } + invisible(list(data_dir = data_dir, reference_dir = reference_dir)) +} + +.hcr_as_df <- function(x) { + if (is.null(x)) { + return(data.frame()) + } + if (inherits(x, "DataFrame")) { + return(as.data.frame(x)) + } + as.data.frame(x, stringsAsFactors = FALSE, check.names = FALSE) +} + +.hcr_write_csv <- function(x, path, row.names = FALSE) { + dir.create(dirname(path), recursive = TRUE, showWarnings = FALSE) + utils::write.csv(x, file = path, row.names = row.names, quote = TRUE, na = "") + normalizePath(path, winslash = "/", mustWork = TRUE) +} + +.hcr_write_json <- function(x, path) { + if (!requireNamespace("jsonlite", quietly = TRUE)) { + stop("Package `jsonlite` is required to write the regression manifest.", call. = FALSE) + } + dir.create(dirname(path), recursive = TRUE, showWarnings = FALSE) + jsonlite::write_json(x, path = path, auto_unbox = TRUE, pretty = TRUE, null = "null") + normalizePath(path, winslash = "/", mustWork = TRUE) +} + +.hcr_sort_df <- function(x, preferred = character()) { + x <- .hcr_as_df(x) + if (nrow(x) == 0) { + return(x) + } + cols <- unique(c(preferred[preferred %in% names(x)], names(x))) + ord_data <- x[cols] + ord <- do.call(order, lapply(ord_data, function(col) { + if (is.numeric(col)) col else tolower(as.character(col)) + })) + x[ord, , drop = FALSE] +} + +.hcr_table_metrics <- function(hc) { + layer_cfg <- .hcr_as_df(hc@config@layer) + experiments <- MultiAssayExperiment::experiments(hc@mae) + layer_rows <- lapply(seq_along(experiments), function(i) { + se <- experiments[[i]] + data.frame( + metric = paste0("layer_", names(experiments)[[i]], "_dimensions"), + value = paste0(nrow(se), " genes x ", ncol(se), " samples"), + stringsAsFactors = FALSE + ) + }) + cfg_rows <- if (nrow(layer_cfg) > 0) { + data.frame( + metric = paste0("configured_layer_", seq_len(nrow(layer_cfg))), + value = apply(layer_cfg, 1, function(row) paste(names(row), row, sep = "=", collapse = "; ")), + stringsAsFactors = FALSE + ) + } else { + data.frame() + } + do.call(rbind, c(layer_rows, list(cfg_rows))) +} + +.hcr_parse_module_genes <- function(gene_n) { + parser <- tryCatch( + get(".hc_parse_genes_from_gene_n", envir = asNamespace("hcocena"), inherits = FALSE), + error = function(e) NULL + ) + if (is.function(parser)) { + out <- tryCatch(parser(gene_n), error = function(e) character()) + return(unique(as.character(out[!is.na(out) & nzchar(out)]))) + } + + gene_n <- as.character(gene_n) + gene_n <- gene_n[!is.na(gene_n)] + if (length(gene_n) == 0) { + return(character()) + } + out <- unlist(strsplit(gene_n, "[#,;|[:space:]]+"), use.names = FALSE) + unique(out[!is.na(out) & nzchar(out)]) +} + +.hcr_choose_split_module <- function(hc) { + cluster <- as.list(hc@integration@cluster) + cluster_info <- .hcr_as_df(cluster[["cluster_information"]]) + label_map <- cluster[["module_label_map"]] + if (nrow(cluster_info) == 0 || !all(c("color", "gene_n") %in% names(cluster_info))) { + stop("Cannot choose a split module: cluster information is missing `color` or `gene_n`.", call. = FALSE) + } + + included <- if ("cluster_included" %in% names(cluster_info)) { + as.character(cluster_info$cluster_included) == "yes" + } else { + rep(TRUE, nrow(cluster_info)) + } + candidates <- cluster_info[included, , drop = FALSE] + if (nrow(candidates) == 0) { + stop("Cannot choose a split module: no included modules found.", call. = FALSE) + } + + gene_counts <- vapply(candidates$gene_n, function(x) length(.hcr_parse_module_genes(x)), integer(1)) + target_idx <- which.max(gene_counts) + target_color <- as.character(candidates$color[[target_idx]]) + target_label <- if (length(label_map) > 0 && target_color %in% names(label_map)) { + as.character(label_map[[target_color]]) + } else { + target_color + } + + list( + module_color = target_color, + module_label = target_label, + gene_count = as.integer(gene_counts[[target_idx]]) + ) +} + +.hcr_empty_df <- function(cols) { + stats::setNames(as.data.frame(rep(list(character()), length(cols))), cols) +} + +.hcr_enrichment_outputs <- function(satellite) { + enrich <- satellite[["enrichments"]] + if (is.null(enrich) || !is.list(enrich)) { + enrich <- list() + } + all <- .hcr_sort_df(enrich[["all_enrichments_all_dbs"]], c("database", "module", "module_label", "ID", "Description", "p.adjust")) + selected <- .hcr_sort_df(enrich[["selected_enrichments_all_dbs"]], c("database", "module", "module_label", "ID", "Description", "p.adjust")) + significant <- .hcr_sort_df(enrich[["significant_enrichments_all_dbs"]], c("database", "module", "module_label", "ID", "Description", "p.adjust")) + if (ncol(all) == 0) all <- .hcr_empty_df(c("database", "module", "module_label", "ID", "Description", "p.adjust")) + if (ncol(selected) == 0) selected <- .hcr_empty_df(c("database", "module", "module_label", "ID", "Description", "p.adjust")) + if (ncol(significant) == 0) significant <- .hcr_empty_df(c("database", "module", "module_label", "ID", "Description", "p.adjust")) + list(all = all, selected = selected, significant = significant) +} + +.hcr_split_outputs <- function(satellite) { + resolution <- satellite[["module_split_resolution_test_last"]] + last <- satellite[["module_split_last"]] + + by_resolution <- if (!is.null(resolution) && is.list(resolution)) { + .hcr_sort_df(resolution[["by_resolution"]], c("module_label", "resolution")) + } else { + data.frame() + } + by_module <- if (!is.null(resolution) && is.list(resolution)) { + .hcr_sort_df(resolution[["by_module"]], c("module_label", "resolution")) + } else { + data.frame() + } + resolved <- if (!is.null(last) && is.list(last)) { + .hcr_sort_df(last[["resolved_modules"]], c("resolved_index", "resolved_label", "resolved_color")) + } else { + data.frame() + } + summary <- if (!is.null(last) && is.list(last)) { + .hcr_sort_df(last[["split_summary"]], c("parent_label", "status")) + } else { + data.frame() + } + + if (ncol(by_resolution) == 0) by_resolution <- .hcr_empty_df(c("module_label", "resolution", "n_submodules")) + if (ncol(by_module) == 0) by_module <- .hcr_empty_df(c("module_label", "resolution", "n_submodules")) + if (ncol(resolved) == 0) resolved <- .hcr_empty_df(c("input", "resolved_index", "resolved_color", "resolved_label", "status")) + if (ncol(summary) == 0) summary <- .hcr_empty_df(c("parent_color", "parent_label", "parent_genes", "n_submodules", "status")) + + list( + by_resolution = by_resolution, + by_module = by_module, + resolved = resolved, + summary = summary + ) +} + +.hcr_visual_text_page <- function(title, lines = character(), footer = NULL) { + grid::grid.newpage() + grid::grid.text( + title, + x = grid::unit(0.05, "npc"), + y = grid::unit(0.94, "npc"), + just = c("left", "top"), + gp = grid::gpar(fontsize = 18, fontface = "bold") + ) + lines <- unlist(lapply(as.character(lines), strwrap, width = 120), use.names = FALSE) + if (length(lines) > 0) { + y <- 0.86 + for (line in lines[seq_len(min(length(lines), 36L))]) { + grid::grid.text( + line, + x = grid::unit(0.05, "npc"), + y = grid::unit(y, "npc"), + just = c("left", "top"), + gp = grid::gpar(fontsize = 9) + ) + y <- y - 0.023 + } + } + if (!is.null(footer)) { + grid::grid.text( + as.character(footer), + x = grid::unit(0.05, "npc"), + y = grid::unit(0.04, "npc"), + just = c("left", "bottom"), + gp = grid::gpar(fontsize = 8, col = "#555555") + ) + } + invisible(NULL) +} + +.hcr_visual_table_page <- function(title, df, n = 18L, footer = NULL) { + df <- .hcr_as_df(df) + if (nrow(df) == 0) { + return(.hcr_visual_text_page(title, "No rows available.", footer = footer)) + } + lines <- utils::capture.output(print(utils::head(df, n), row.names = FALSE, right = FALSE)) + .hcr_visual_text_page(title, lines, footer = footer) +} + +.hcr_visual_parameters_df <- function(parameters) { + if (is.null(parameters) || length(parameters) == 0) { + return(data.frame()) + } + values <- vapply(parameters, .hcr_value_label, FUN.VALUE = character(1), USE.NAMES = FALSE) + data.frame( + parameter = names(parameters), + value = values, + stringsAsFactors = FALSE + ) +} + +.hcr_value_label <- function(value) { + if (is.null(value)) { + return("NULL") + } + if (is.atomic(value) && length(value) <= 12) { + return(paste(as.character(value), collapse = ", ")) + } + paste(utils::capture.output(str(value, give.attr = FALSE)), collapse = " ") +} + +.hcr_values_equal <- function(value, default) { + if (is.null(value) && is.null(default)) { + return(TRUE) + } + if (is.null(value) || is.null(default)) { + return(FALSE) + } + if (is.numeric(value) && is.numeric(default) && length(value) == length(default)) { + return(isTRUE(all.equal(as.numeric(value), as.numeric(default), tolerance = 1e-12))) + } + identical(as.character(value), as.character(default)) +} + +.hcr_parameter_diffs_df <- function(parameters, defaults = list(), scope = "") { + if (is.null(parameters) || length(parameters) == 0) { + return(data.frame()) + } + rows <- lapply(names(parameters), function(parameter) { + value <- parameters[[parameter]] + default <- defaults[[parameter]] + if (.hcr_values_equal(value, default)) { + return(NULL) + } + data.frame( + scope = scope, + parameter = parameter, + value = .hcr_value_label(value), + default = .hcr_value_label(default), + stringsAsFactors = FALSE + ) + }) + rows <- rows[!vapply(rows, is.null, logical(1))] + if (length(rows) == 0) { + return(data.frame()) + } + do.call(rbind, rows) +} + +.hcr_heatmap_plot_defaults <- function() { + list( + module_label_preset = "balanced", + module_label_fontsize = NULL, + module_label_pt_size = NULL, + module_box_width_cm = NULL, + gene_count_mode = "text" + ) +} + +.hcr_enrichment_plot_defaults <- function() { + list( + gene_sets = "Hallmark", + top = 5, + consistent_terms = TRUE, + cluster_columns = FALSE, + store_panel_objects = "auto", + heatmap_module_label_fontsize = NULL + ) +} + +.hcr_plot_catalog_df <- function(plot_variants) { + if (is.null(plot_variants) || length(plot_variants) == 0) { + return(data.frame()) + } + rows <- lapply(plot_variants, function(variant) { + data.frame( + plot = .hcr_scalar(variant$plot, ""), + variant = .hcr_scalar(variant$variant, ""), + pdf_file = .hcr_scalar(variant$pdf_file, ""), + png_file = .hcr_scalar(variant$png_file, ""), + stage = .hcr_scalar(variant$stage, ""), + stringsAsFactors = FALSE + ) + }) + do.call(rbind, rows) +} + +.hcr_plot_parameter_diffs_df <- function(plot_variants, font_probe = data.frame()) { + if (is.null(plot_variants) || length(plot_variants) == 0) { + return(data.frame()) + } + rows <- list() + for (variant in plot_variants) { + variant_name <- paste( + c(.hcr_scalar(variant$plot, ""), .hcr_scalar(variant$variant, "")), + collapse = " | " + ) + params <- variant$parameters + defaults <- variant$defaults + diff <- .hcr_parameter_diffs_df(params, defaults, scope = variant_name) + if (nrow(diff) > 0) { + diff$source <- "non-default parameter" + rows[[length(rows) + 1L]] <- diff + } + } + + font_probe_scope <- "Cluster heatmap | module-label font-size probe" + font_probe_idx <- which(vapply(plot_variants, function(variant) { + identical(.hcr_scalar(variant$plot, ""), "Cluster heatmap") && + grepl("font-size probe", .hcr_scalar(variant$variant, ""), ignore.case = TRUE) + }, FUN.VALUE = logical(1))) + if (length(font_probe_idx) > 0) { + font_probe_variant <- plot_variants[[font_probe_idx[[1]]]] + font_probe_scope <- paste( + c(.hcr_scalar(font_probe_variant$plot, ""), .hcr_scalar(font_probe_variant$variant, "")), + collapse = " | " + ) + } + + if (nrow(font_probe) > 0 && + "module_label_auto_fit_shrunk" %in% names(font_probe) && + isTRUE(tolower(as.character(font_probe$module_label_auto_fit_shrunk[[1]])) %in% c("true", "1", "yes"))) { + reason <- paste( + c( + if ("module_label_auto_fit_width_limited" %in% names(font_probe) && + tolower(as.character(font_probe$module_label_auto_fit_width_limited[[1]])) %in% c("true", "1", "yes")) "width", + if ("module_label_auto_fit_height_limited" %in% names(font_probe) && + tolower(as.character(font_probe$module_label_auto_fit_height_limited[[1]])) %in% c("true", "1", "yes")) "height" + ), + collapse = "+" + ) + if (!nzchar(reason)) { + reason <- "fit guard" + } + rows[[length(rows) + 1L]] <- data.frame( + scope = font_probe_scope, + parameter = "effective_module_label_font_size_pt", + value = paste0( + .hcr_scalar(font_probe$effective_module_label_pt_size_min_pt, ""), + " to ", + .hcr_scalar(font_probe$effective_module_label_pt_size_max_pt, "") + ), + default = paste0("requested max ", .hcr_scalar(font_probe$requested_module_label_pt_size_max_pt, "")), + source = paste0("runtime auto-fit (", reason, ")"), + stringsAsFactors = FALSE + ) + } + + if (length(rows) == 0) { + return(data.frame()) + } + do.call(rbind, rows) +} + +.hcr_run_parameter_defaults <- function(mode) { + c( + list( + mode = mode, + seed = 168575, + compare_reference = TRUE, + update_reference = FALSE, + tolerance = 1e-8 + ), + .hcr_mode_config(mode) + ) +} + +.hcr_plot_file <- function(path, base_dir) { + rel <- normalizePath(path, winslash = "/", mustWork = FALSE) + root <- normalizePath(base_dir, winslash = "/", mustWork = FALSE) + prefix <- paste0(root, "/") + if (startsWith(rel, prefix)) { + rel <- substring(rel, nchar(prefix) + 1L) + } + rel +} + +.hcr_plot_variant <- function(plot, + variant, + pdf_file, + base_dir, + stage, + parameters, + defaults) { + png_file <- sub("\\.pdf$", ".png", pdf_file, ignore.case = TRUE) + list( + plot = plot, + variant = variant, + pdf_file = .hcr_plot_file(pdf_file, base_dir), + png_file = if (file.exists(png_file)) .hcr_plot_file(png_file, base_dir) else "", + stage = stage, + parameters = parameters, + defaults = defaults + ) +} + +.hcr_run_parameter_diffs_df <- function(parameters, mode) { + keep <- setdiff(names(parameters), c("data_dir", "output_dir", "reference_dir", "package_source", "package_version")) + .hcr_parameter_diffs_df(parameters[keep], .hcr_run_parameter_defaults(mode), scope = "runner") +} + +.hcr_visual_parameter_table <- function(title, df, empty_message = "No non-default parameters.") { + df <- .hcr_as_df(df) + if (nrow(df) == 0) { + return(.hcr_visual_text_page(title, empty_message)) + } + .hcr_visual_table_page(title, df, n = 48L) +} + +.hcr_is_abs_path <- function(path) { + path <- as.character(path) + grepl("^[A-Za-z]:[/\\\\]", path) || grepl("^[/\\\\]{2}", path) || startsWith(path, "/") +} + +.hcr_abs_file <- function(path, base_dir) { + path <- .hcr_scalar(path, "") + if (!nzchar(path)) { + return("") + } + if (.hcr_is_abs_path(path)) { + return(normalizePath(path, winslash = "/", mustWork = FALSE)) + } + normalizePath(file.path(base_dir, path), winslash = "/", mustWork = FALSE) +} + +.hcr_plot_png_file <- function(variant, output_dir) { + png_file <- .hcr_scalar(variant$png_file, "") + if (!nzchar(png_file)) { + pdf_file <- .hcr_scalar(variant$pdf_file, "") + if (nzchar(pdf_file)) { + png_file <- sub("\\.pdf$", ".png", pdf_file, ignore.case = TRUE) + } + } + .hcr_abs_file(png_file, output_dir) +} + +.hcr_plot_detail_lines <- function(variant, output_dir, font_probe = data.frame()) { + pdf_file <- .hcr_scalar(variant$pdf_file, "") + png_file <- .hcr_scalar(variant$png_file, "") + if (!nzchar(png_file) && nzchar(pdf_file)) { + png_file <- sub("\\.pdf$", ".png", pdf_file, ignore.case = TRUE) + } + diff <- .hcr_parameter_diffs_df(variant$parameters, variant$defaults) + parameter_line <- if (nrow(diff) == 0) { + "Non-default parameters: none" + } else { + paste0( + "Non-default parameters: ", + paste0(diff$parameter, "=", diff$value, " (default ", diff$default, ")", collapse = "; ") + ) + } + lines <- c( + paste0("Stage: ", .hcr_scalar(variant$stage, "")), + paste0("PDF: ", pdf_file), + paste0("PNG: ", png_file), + parameter_line + ) + + is_font_probe <- grepl("font-size probe", .hcr_scalar(variant$variant, ""), ignore.case = TRUE) + if (is_font_probe && nrow(font_probe) > 0) { + lines <- c( + lines, + paste0( + "Runtime module label size: ", + .hcr_scalar(font_probe$effective_module_label_pt_size_min_pt, ""), + " to ", + .hcr_scalar(font_probe$effective_module_label_pt_size_max_pt, ""), + " pt (requested max ", + .hcr_scalar(font_probe$requested_module_label_pt_size_max_pt, ""), + " pt; auto-fit shrunk=", + .hcr_scalar(font_probe$module_label_auto_fit_shrunk, ""), + ")" + ) + ) + } + lines +} + +.hcr_visual_plot_page <- function(title, png_file, lines = character(), footer = NULL) { + if (!nzchar(.hcr_scalar(png_file, "")) || !file.exists(png_file)) { + return(.hcr_visual_text_page( + title, + c("Plot image is unavailable.", paste0("PNG: ", png_file), lines), + footer = footer + )) + } + if (!requireNamespace("png", quietly = TRUE)) { + return(.hcr_visual_text_page( + title, + c("Package `png` is not installed, so the plot image cannot be embedded.", paste0("PNG: ", png_file), lines), + footer = footer + )) + } + + image <- tryCatch( + png::readPNG(png_file), + error = function(e) e + ) + if (inherits(image, "error")) { + return(.hcr_visual_text_page( + title, + c("Plot image could not be read.", paste0("PNG: ", png_file), paste0("Error: ", conditionMessage(image)), lines), + footer = footer + )) + } + + grid::grid.newpage() + grid::grid.text( + title, + x = grid::unit(0.05, "npc"), + y = grid::unit(0.985, "npc"), + just = c("left", "top"), + gp = grid::gpar(fontsize = 11, fontface = "bold") + ) + + detail_lines <- unlist(lapply(as.character(lines), strwrap, width = 150), use.names = FALSE) + detail_lines <- detail_lines[seq_len(min(length(detail_lines), 6L))] + + detail_bottom <- if (is.null(footer)) 0.035 else 0.055 + detail_line_height <- 0.017 + detail_height <- if (length(detail_lines) > 0) { + (length(detail_lines) * detail_line_height) + 0.014 + } else { + 0 + } + image_top <- 0.955 + image_bottom <- detail_bottom + detail_height + max_w <- 0.96 + max_h <- max(0.25, image_top - image_bottom) + aspect <- dim(image)[[2]] / dim(image)[[1]] + if (aspect >= max_w / max_h) { + width <- max_w + height <- max_w / aspect + } else { + height <- max_h + width <- max_h * aspect + } + grid::grid.raster( + image, + x = grid::unit(0.5, "npc"), + y = grid::unit(image_bottom + (max_h / 2), "npc"), + width = grid::unit(width, "npc"), + height = grid::unit(height, "npc"), + interpolate = TRUE + ) + + if (length(detail_lines) > 0) { + y <- detail_bottom + (length(detail_lines) * detail_line_height) + for (line in detail_lines) { + grid::grid.text( + line, + x = grid::unit(0.05, "npc"), + y = grid::unit(y, "npc"), + just = c("left", "top"), + gp = grid::gpar(fontsize = 6.4, col = "#333333") + ) + y <- y - detail_line_height + } + } + + if (!is.null(footer)) { + grid::grid.text( + as.character(footer), + x = grid::unit(0.05, "npc"), + y = grid::unit(0.028, "npc"), + just = c("left", "bottom"), + gp = grid::gpar(fontsize = 6.8, col = "#555555") + ) + } + invisible(NULL) +} + +.hcr_visual_plot_pages <- function(plot_variants, output_dir, font_probe = data.frame()) { + if (is.null(plot_variants) || length(plot_variants) == 0) { + return(invisible(NULL)) + } + footer <- normalizePath(output_dir, winslash = "/", mustWork = FALSE) + for (variant in plot_variants) { + title <- paste( + c("Plot image", .hcr_scalar(variant$plot, ""), .hcr_scalar(variant$variant, "")), + collapse = " - " + ) + .hcr_visual_plot_page( + title = title, + png_file = .hcr_plot_png_file(variant, output_dir), + lines = .hcr_plot_detail_lines(variant, output_dir, font_probe = font_probe), + footer = footer + ) + } + invisible(NULL) +} + +.hcr_open_visual_report_pdf <- function(report_file, width = 13, height = 9.5) { + candidates <- report_file + fallback <- file.path( + dirname(report_file), + paste0( + tools::file_path_sans_ext(basename(report_file)), + "_", + format(Sys.time(), "%Y%m%d_%H%M%S"), + ".pdf" + ) + ) + candidates <- unique(c(candidates, fallback)) + errors <- character() + for (candidate in candidates) { + ok <- tryCatch({ + grDevices::pdf(candidate, width = width, height = height, onefile = TRUE) + TRUE + }, error = function(e) { + errors <<- c(errors, paste0(candidate, ": ", conditionMessage(e))) + FALSE + }) + if (isTRUE(ok)) { + return(normalizePath(candidate, winslash = "/", mustWork = FALSE)) + } + } + stop( + "Could not open a visual-check PDF device. Tried:\n", + paste(errors, collapse = "\n"), + call. = FALSE + ) +} + +.hcr_visual_parameters_df_legacy <- function(parameters) { + if (is.null(parameters) || length(parameters) == 0) { + return(data.frame()) + } + values <- vapply(parameters, function(value) { + if (is.null(value)) { + return("NULL") + } + if (is.atomic(value) && length(value) <= 12) { + return(paste(as.character(value), collapse = ", ")) + } + paste(utils::capture.output(str(value, give.attr = FALSE)), collapse = " ") + }, FUN.VALUE = character(1), USE.NAMES = FALSE) + data.frame( + parameter = names(parameters), + value = values, + stringsAsFactors = FALSE + ) +} + +.hcr_write_visual_report <- function(hc, + output_dir, + mode, + font_probe_file, + split_target = NULL, + parameters = NULL, + plot_variants = NULL) { + report_file <- file.path(output_dir, paste0("visual_check_report_", mode, ".pdf")) + dir.create(dirname(report_file), recursive = TRUE, showWarnings = FALSE) + + read_export <- function(name) { + path <- file.path(output_dir, name) + if (!file.exists(path)) { + return(data.frame()) + } + utils::read.csv(path, check.names = FALSE, stringsAsFactors = FALSE) + } + + qc <- read_export("qc_summary.csv") + split_summary <- read_export("split_summary.csv") + split_resolution <- read_export("split_resolution_summary.csv") + enrichment_selected <- read_export("enrichment_selected.csv") + font_probe <- read_export("font_probe.csv") + plot_catalog <- .hcr_plot_catalog_df(plot_variants) + plot_parameter_diffs <- .hcr_plot_parameter_diffs_df(plot_variants, font_probe = font_probe) + run_parameter_diffs <- .hcr_run_parameter_diffs_df(parameters, mode = mode) + + report_file <- .hcr_open_visual_report_pdf(report_file, width = 13, height = 9.5) + on.exit(grDevices::dev.off(), add = TRUE) + + qc_lines <- c( + paste0("Mode: ", mode), + paste0("Generated: ", format(Sys.time(), "%Y-%m-%d %H:%M:%S")), + if (!is.null(split_target)) { + paste0( + "Split target: ", split_target$module_label, + " (", split_target$gene_count, " genes)" + ) + } else { + "Split target: none" + }, + "", + utils::capture.output(print(qc, row.names = FALSE, right = FALSE)) + ) + .hcr_visual_text_page( + title = "hCoCena real-data visual check", + lines = qc_lines, + footer = normalizePath(output_dir, winslash = "/", mustWork = FALSE) + ) + + .hcr_visual_table_page("Plot variants", plot_catalog, n = 32L) + + .hcr_visual_parameter_table( + "Plot parameter deviations", + plot_parameter_diffs, + empty_message = "No plot parameters differ from their defaults." + ) + + .hcr_visual_parameter_table( + "Run parameter deviations", + run_parameter_diffs, + empty_message = "No runner parameters differ from the mode defaults." + ) + + .hcr_visual_plot_pages(plot_variants, output_dir, font_probe = font_probe) + + .hcr_visual_table_page("Module split summary", split_summary) + + if (requireNamespace("ggplot2", quietly = TRUE) && nrow(split_resolution) > 0 && + all(c("resolution", "total_created_submodules", "resulting_included_modules") %in% names(split_resolution))) { + p <- ggplot2::ggplot(split_resolution, ggplot2::aes(x = resolution)) + + ggplot2::geom_line(ggplot2::aes(y = total_created_submodules, color = "created submodules"), linewidth = 0.8) + + ggplot2::geom_point(ggplot2::aes(y = total_created_submodules, color = "created submodules"), size = 2) + + ggplot2::geom_line(ggplot2::aes(y = resulting_included_modules, color = "included modules"), linewidth = 0.8) + + ggplot2::geom_point(ggplot2::aes(y = resulting_included_modules, color = "included modules"), size = 2) + + ggplot2::labs( + title = "Split resolution probe", + x = "Resolution", + y = "Count", + color = NULL + ) + + ggplot2::theme_minimal(base_size = 12) + print(p) + } else { + .hcr_visual_table_page("Split resolution probe", split_resolution) + } + + if (requireNamespace("ggplot2", quietly = TRUE) && nrow(enrichment_selected) > 0 && + all(c("database", "module_label", "term", "p.adjust") %in% names(enrichment_selected))) { + enrich_plot <- enrichment_selected + enrich_plot$p_adjust_num <- suppressWarnings(as.numeric(enrich_plot$p.adjust)) + enrich_plot <- enrich_plot[is.finite(enrich_plot$p_adjust_num), , drop = FALSE] + enrich_plot <- enrich_plot[order(enrich_plot$p_adjust_num), , drop = FALSE] + enrich_plot <- utils::head(enrich_plot, 18L) + enrich_plot$label <- paste0(enrich_plot$module_label, " | ", enrich_plot$database, " | ", enrich_plot$term) + enrich_plot$label <- factor(enrich_plot$label, levels = rev(enrich_plot$label)) + p <- ggplot2::ggplot(enrich_plot, ggplot2::aes(x = -log10(p_adjust_num), y = label)) + + ggplot2::geom_col(fill = "#4C78A8") + + ggplot2::labs( + title = "Selected enrichment terms", + x = "-log10(adjusted p)", + y = NULL + ) + + ggplot2::theme_minimal(base_size = 9) + + ggplot2::theme(plot.title = ggplot2::element_text(size = 14, face = "bold")) + print(p) + } else { + .hcr_visual_table_page("Selected enrichment terms", enrichment_selected) + } + + .hcr_visual_table_page("Module label font-size probe", font_probe, n = 8L) + + cluster <- as.list(hc@integration@cluster) + heatmap_obj <- if (!is.null(cluster[["heatmap_cluster"]])) { + cluster[["heatmap_cluster"]] + } else { + cluster[["heatmap_cluster_raw"]] + } + if (!is.null(heatmap_obj) && requireNamespace("ComplexHeatmap", quietly = TRUE)) { + ok <- tryCatch({ + ComplexHeatmap::draw(heatmap_obj, newpage = TRUE, merge_legends = TRUE) + TRUE + }, error = function(e) { + .hcr_visual_text_page( + "Cluster heatmap", + c( + "The heatmap object could not be redrawn in the visual report.", + paste0("Standalone heatmap file: ", file.path("hcocena", basename(font_probe_file))), + paste0("Error: ", conditionMessage(e)) + ) + ) + FALSE + }) + invisible(ok) + } else { + .hcr_visual_text_page( + "Cluster heatmap", + c( + "No drawable heatmap object was stored.", + paste0("Standalone heatmap file: ", file.path("hcocena", basename(font_probe_file))) + ) + ) + } + + normalizePath(report_file, winslash = "/", mustWork = TRUE) +} + +.hcr_export_outputs <- function(hc, output_dir, mode, font_probe_file, plot_variants = NULL) { + cluster <- as.list(hc@integration@cluster) + satellite <- as.list(hc@satellite) + enrichment <- .hcr_enrichment_outputs(satellite) + split <- .hcr_split_outputs(satellite) + cluster_info <- .hcr_sort_df(cluster[["cluster_information"]], c("cluster", "color", "gene_n")) + edge_list <- .hcr_sort_df(hc@integration@combined_edgelist, c("V1", "V2", "weight")) + gfc <- .hcr_sort_df(hc@integration@gfc, c("Gene")) + module_gene_src <- if (!is.null(satellite[["module_gene_list"]])) { + satellite[["module_gene_list"]] + } else { + cluster[["module_gene_list"]] + } + module_gene_list <- .hcr_sort_df(module_gene_src, c("module", "module_color", "gene")) + if (ncol(module_gene_list) == 0) { + module_gene_list <- data.frame(module = character(), module_color = character(), gene = character()) + } + + label_map <- cluster[["module_label_map"]] + label_map_df <- if (length(label_map) > 0) { + data.frame( + module_color = names(label_map), + module_label = as.character(label_map), + stringsAsFactors = FALSE + ) + } else { + data.frame(module_color = character(), module_label = character()) + } + label_map_df <- .hcr_sort_df(label_map_df, c("module_label", "module_color")) + + heatmap_matrix <- cluster[["heatmap_matrix"]] + heatmap_df <- if (!is.null(heatmap_matrix) && length(heatmap_matrix) > 0) { + data.frame(module = rownames(heatmap_matrix), as.data.frame(heatmap_matrix, check.names = FALSE), check.names = FALSE) + } else { + data.frame() + } + + graph <- hc@integration@graph + graph_metrics <- data.frame( + metric = c("nodes", "edges"), + value = c( + if (inherits(graph, "igraph")) igraph::vcount(graph) else NA_integer_, + if (inherits(graph, "igraph")) igraph::ecount(graph) else NA_integer_ + ), + stringsAsFactors = FALSE + ) + module_metrics <- data.frame( + metric = c( + "included_modules", + "cluster_rows", + "module_label_fontsize", + "module_label_pt_size", + "module_label_pt_size_effective_pt_min", + "module_label_pt_size_effective_pt_max", + "module_label_pt_size_requested_pt_max", + "module_label_auto_fit_shrunk", + "module_label_auto_fit_width_limited", + "module_label_auto_fit_height_limited", + "module_box_width_cm", + "split_summary_rows", + "enrichment_all_rows", + "enrichment_selected_rows", + "enrichment_significant_rows", + "heatmap_rows", + "heatmap_columns" + ), + value = c( + if (nrow(cluster_info) > 0 && "cluster_included" %in% names(cluster_info)) { + length(unique(cluster_info$color[cluster_info$cluster_included == "yes"])) + } else { + NA_integer_ + }, + nrow(cluster_info), + .hcr_scalar(cluster[["module_label_fontsize"]], NA_character_), + .hcr_scalar(cluster[["module_label_pt_size"]], NA_character_), + .hcr_scalar(cluster[["module_label_pt_size_effective_pt_min"]], NA_character_), + .hcr_scalar(cluster[["module_label_pt_size_effective_pt_max"]], NA_character_), + .hcr_scalar(cluster[["module_label_pt_size_requested_pt_max"]], NA_character_), + .hcr_scalar(cluster[["module_label_auto_fit_shrunk"]], NA_character_), + .hcr_scalar(cluster[["module_label_auto_fit_width_limited"]], NA_character_), + .hcr_scalar(cluster[["module_label_auto_fit_height_limited"]], NA_character_), + .hcr_scalar(cluster[["module_box_width_cm"]], NA_character_), + nrow(split$summary), + nrow(enrichment$all), + nrow(enrichment$selected), + nrow(enrichment$significant), + nrow(heatmap_df), + max(0L, ncol(heatmap_df) - 1L) + ), + stringsAsFactors = FALSE + ) + font_probe <- data.frame( + requested_module_label_fontsize = .hcr_scalar(cluster[["module_label_fontsize"]], NA_character_), + requested_module_label_pt_size = .hcr_scalar(cluster[["module_label_pt_size"]], NA_character_), + requested_module_box_width_cm = .hcr_scalar(cluster[["module_box_width_cm"]], NA_character_), + requested_module_label_pt_size_max_pt = .hcr_scalar(cluster[["module_label_pt_size_requested_pt_max"]], NA_character_), + effective_module_label_pt_size_min_pt = .hcr_scalar(cluster[["module_label_pt_size_effective_pt_min"]], NA_character_), + effective_module_label_pt_size_max_pt = .hcr_scalar(cluster[["module_label_pt_size_effective_pt_max"]], NA_character_), + module_label_auto_fit_shrunk = .hcr_scalar(cluster[["module_label_auto_fit_shrunk"]], NA_character_), + module_label_auto_fit_width_limited = .hcr_scalar(cluster[["module_label_auto_fit_width_limited"]], NA_character_), + module_label_auto_fit_height_limited = .hcr_scalar(cluster[["module_label_auto_fit_height_limited"]], NA_character_), + effective_size_control = "module_label_pt_size_with_box_fit_guard", + heatmap_file = file.path("hcocena", basename(font_probe_file)), + heatmap_file_exists = file.exists(font_probe_file), + stringsAsFactors = FALSE + ) + plot_catalog <- .hcr_plot_catalog_df(plot_variants) + plot_parameter_diffs <- .hcr_plot_parameter_diffs_df(plot_variants, font_probe = font_probe) + qc <- rbind(.hcr_table_metrics(hc), graph_metrics, module_metrics) + + files <- list( + qc_summary = .hcr_write_csv(qc, file.path(output_dir, "qc_summary.csv")), + combined_edgelist = .hcr_write_csv(edge_list, file.path(output_dir, "combined_edgelist.csv")), + gfc_matrix = .hcr_write_csv(gfc, file.path(output_dir, "gfc_matrix.csv")), + cluster_information = .hcr_write_csv(cluster_info, file.path(output_dir, "cluster_information.csv")), + module_label_map = .hcr_write_csv(label_map_df, file.path(output_dir, "module_label_map.csv")), + module_gene_list = .hcr_write_csv(module_gene_list, file.path(output_dir, "module_gene_list.csv")), + heatmap_matrix = .hcr_write_csv(heatmap_df, file.path(output_dir, "heatmap_matrix.csv")), + font_probe = .hcr_write_csv(font_probe, file.path(output_dir, "font_probe.csv")), + plot_variants = .hcr_write_csv(plot_catalog, file.path(output_dir, "plot_variants.csv")), + plot_variant_parameters = .hcr_write_csv(plot_parameter_diffs, file.path(output_dir, "plot_variant_parameters.csv")), + split_resolution_summary = .hcr_write_csv(split$by_resolution, file.path(output_dir, "split_resolution_summary.csv")), + split_resolution_by_module = .hcr_write_csv(split$by_module, file.path(output_dir, "split_resolution_by_module.csv")), + split_resolved_modules = .hcr_write_csv(split$resolved, file.path(output_dir, "split_resolved_modules.csv")), + split_summary = .hcr_write_csv(split$summary, file.path(output_dir, "split_summary.csv")), + enrichment_all = .hcr_write_csv(enrichment$all, file.path(output_dir, "enrichment_all.csv")), + enrichment_selected = .hcr_write_csv(enrichment$selected, file.path(output_dir, "enrichment_selected.csv")), + enrichment_significant = .hcr_write_csv(enrichment$significant, file.path(output_dir, "enrichment_significant.csv")) + ) + + list( + files = files, + metrics = qc, + font_probe = font_probe, + plot_variants = plot_catalog, + plot_variant_parameters = plot_parameter_diffs + ) +} + +.hcr_compare_table <- function(current_path, reference_path, tolerance = 1e-8) { + if (!file.exists(reference_path)) { + return(list(status = "missing_reference", details = "reference file not found")) + } + cur <- utils::read.csv(current_path, check.names = FALSE, stringsAsFactors = FALSE) + ref <- utils::read.csv(reference_path, check.names = FALSE, stringsAsFactors = FALSE) + if (!identical(dim(cur), dim(ref))) { + return(list(status = "different", details = paste0("dimension ", paste(dim(cur), collapse = "x"), " != ", paste(dim(ref), collapse = "x")))) + } + if (!identical(names(cur), names(ref))) { + return(list(status = "different", details = "column names differ")) + } + for (nm in names(cur)) { + cx <- cur[[nm]] + rx <- ref[[nm]] + c_num <- suppressWarnings(as.numeric(cx)) + r_num <- suppressWarnings(as.numeric(rx)) + numeric_like <- !all(is.na(c_num)) || !all(is.na(r_num)) + if (numeric_like) { + bad <- xor(is.na(c_num), is.na(r_num)) | + (!is.na(c_num) & !is.na(r_num) & abs(c_num - r_num) > tolerance) + if (any(bad)) { + return(list(status = "different", details = paste0("numeric column differs: ", nm))) + } + } else if (!identical(as.character(cx), as.character(rx))) { + return(list(status = "different", details = paste0("character column differs: ", nm))) + } + } + list(status = "ok", details = "") +} + +.hcr_compare_outputs <- function(output_dir, reference_dir, tolerance = 1e-8) { + files <- c( + "qc_summary.csv", + "combined_edgelist.csv", + "gfc_matrix.csv", + "cluster_information.csv", + "module_label_map.csv", + "module_gene_list.csv", + "heatmap_matrix.csv", + "font_probe.csv", + "plot_variants.csv", + "plot_variant_parameters.csv", + "split_resolution_summary.csv", + "split_resolution_by_module.csv", + "split_resolved_modules.csv", + "split_summary.csv", + "enrichment_all.csv", + "enrichment_selected.csv", + "enrichment_significant.csv" + ) + rows <- lapply(files, function(file) { + res <- .hcr_compare_table( + current_path = file.path(output_dir, file), + reference_path = file.path(reference_dir, file), + tolerance = tolerance + ) + data.frame(file = file, status = res$status, details = res$details, stringsAsFactors = FALSE) + }) + do.call(rbind, rows) +} + +.hcr_update_reference <- function(output_dir, reference_dir) { + dir.create(reference_dir, recursive = TRUE, showWarnings = FALSE) + files <- list.files(output_dir, pattern = "\\.(csv|json)$", full.names = TRUE) + ok <- file.copy(files, reference_dir, overwrite = TRUE) + if (!all(ok)) { + stop("Could not update all reference files in: ", reference_dir, call. = FALSE) + } + invisible(normalizePath(reference_dir, winslash = "/", mustWork = TRUE)) +} + +run_hcocena_realdata_regression <- function(mode = c("quick", "full"), + data_dir = NULL, + output_dir = NULL, + reference_dir = NULL, + update_reference = FALSE, + compare_reference = TRUE, + load_source = TRUE, + seed = 168575, + tolerance = 1e-8, + repo_root = NULL) { + mode <- match.arg(mode) + repo_root <- if (is.null(repo_root)) .hcr_repo_root() else normalizePath(repo_root, winslash = "/", mustWork = TRUE) + cfg <- .hcr_mode_config(mode) + data_dir <- .hcr_scalar(data_dir, .hcr_default_data_dir(repo_root)) + output_dir <- .hcr_scalar(output_dir, file.path(repo_root, "realdata-output", mode)) + reference_dir <- .hcr_scalar(reference_dir, .hcr_default_reference_dir(repo_root, mode)) + + .hcr_check_inputs(data_dir, reference_dir) + data_dir <- normalizePath(data_dir, winslash = "/", mustWork = TRUE) + output_dir <- normalizePath(output_dir, winslash = "/", mustWork = FALSE) + reference_dir <- normalizePath(reference_dir, winslash = "/", mustWork = FALSE) + + dir.create(output_dir, recursive = TRUE, showWarnings = FALSE) + hc_output_dir <- file.path(output_dir, "hcocena") + dir.create(hc_output_dir, recursive = TRUE, showWarnings = FALSE) + + package_source <- .hcr_load_package(repo_root, load_source = load_source) + set.seed(as.integer(seed)) + + reference_files <- file.path(repo_root, "docker", "reference_files") + if (!dir.exists(reference_files)) { + reference_files <- file.path(repo_root, "inst", "extdata") + } + + message("Running hCoCena real-data regression: mode=", mode) + message("Data dir: ", data_dir) + message("Output dir: ", output_dir) + + hc <- hcocena::hc_init() + hc <- hcocena::hc_set_paths( + hc, + dir_count_data = .hcr_slash(data_dir), + dir_annotation = .hcr_slash(data_dir), + dir_reference_files = .hcr_slash(reference_files), + dir_output = .hcr_slash(hc_output_dir) + ) + hc <- hcocena::hc_define_layers( + hc, + data_sets = list( + RNA_Seq = c("data_seq_processed.txt", "annotation_seq.txt"), + Array = c("data_array_processed.txt", "annotation_array.txt") + ) + ) + hc <- hcocena::hc_read_data( + hc, + sep_counts = "\t", + sep_anno = "\t", + gene_symbol_col = "SYMBOL", + sample_col = "SampleID", + count_has_rn = FALSE, + anno_has_rn = FALSE, + project_folder = "" + ) + + if (isTRUE(cfg$read_supplementary)) { + hc <- hcocena::hc_set_supp_files( + hc, + Tf = "TFcat.txt", + Hallmark = "h.all.v2023.1.Hs.symbols.gmt", + Go = "c5.go.bp.v2023.1.Hs.symbols.gmt", + Kegg = "c2.cp.kegg.v2023.1.Hs.symbols.gmt", + Reactome = "c2.cp.reactome.v2023.1.Hs.symbols.gmt" + ) + hc <- hcocena::hc_read_supplementary(hc) + } + + hc <- hcocena::hc_set_global_settings( + hc, + organism = "human", + control_keyword = "baseline", + variable_of_interest = "merged", + min_nodes_number_for_network = cfg$min_nodes_number_for_network, + min_nodes_number_for_cluster = cfg$min_nodes_number_for_cluster, + range_GFC = 2.0, + layout_algorithm = "layout_with_fr", + data_in_log = TRUE + ) + hc <- hcocena::hc_set_layer_settings( + hc, + top_var = cfg$top_var, + min_corr = cfg$min_corr, + range_cutoff_length = cfg$range_cutoff_length, + print_distribution_plots = c(FALSE, FALSE) + ) + hc <- hcocena::hc_run_expression_analysis_1(hc, corr_method = "pearson") + hc <- hcocena::hc_set_cutoff(hc, cutoff_vector = cfg$cutoffs, verbose = FALSE) + hc <- hcocena::hc_run_expression_analysis_2(hc, plot_HM = isTRUE(cfg$plot_layer_heatmaps)) + hc <- hcocena::hc_build_integrated_network(hc, mode = "u", multi_edges = "min") + hc <- hcocena::hc_cluster_calculation( + hc, + cluster_algo = cfg$cluster_algo, + no_of_iterations = cfg$no_of_iterations, + resolution = cfg$resolution + ) + + plot_variants <- list() + heatmap_plot_defaults <- .hcr_heatmap_plot_defaults() + default_heatmap <- file.path(hc_output_dir, paste0("cluster_heatmap_", mode, ".pdf")) + default_heatmap_parameters <- list( + module_label_preset = "balanced", + module_label_fontsize = NULL, + module_label_pt_size = NULL, + module_box_width_cm = NULL, + gene_count_mode = "text" + ) + hc <- hcocena::hc_plot_cluster_heatmap( + hc, + file_name = basename(default_heatmap), + module_label_preset = default_heatmap_parameters$module_label_preset, + module_label_fontsize = default_heatmap_parameters$module_label_fontsize, + module_label_pt_size = default_heatmap_parameters$module_label_pt_size, + module_box_width_cm = default_heatmap_parameters$module_box_width_cm, + gene_count_mode = default_heatmap_parameters$gene_count_mode, + return_HM = FALSE + ) + plot_variants[[length(plot_variants) + 1L]] <- .hcr_plot_variant( + plot = "Cluster heatmap", + variant = "baseline before module split", + pdf_file = default_heatmap, + base_dir = output_dir, + stage = "post clustering", + parameters = default_heatmap_parameters, + defaults = heatmap_plot_defaults + ) + + split_target <- NULL + if (isTRUE(cfg$run_module_split)) { + split_target <- .hcr_choose_split_module(hc) + message( + "Testing module split on ", + split_target$module_label, + " (", split_target$gene_count, " genes)." + ) + hc <- hcocena::hc_split_modules( + hc, + modules = split_target$module_label, + cluster_algo = cfg$split_cluster_algo, + no_of_iterations = 1, + resolution = cfg$split_resolution, + resolution_grid = cfg$split_resolution_grid, + resolution_test_only = TRUE, + seed = seed, + min_submodule_size = cfg$split_min_submodule_size, + verbose = FALSE + ) + hc <- hcocena::hc_split_modules( + hc, + modules = split_target$module_label, + cluster_algo = cfg$split_cluster_algo, + no_of_iterations = 1, + resolution = cfg$split_resolution, + resolution_grid = cfg$split_resolution_grid, + resolution_test_only = FALSE, + seed = seed, + min_submodule_size = cfg$split_min_submodule_size, + verbose = FALSE + ) + } + + font_probe_suffix <- if (isTRUE(cfg$run_module_split)) { + "_post_split_module_label_size_probe.pdf" + } else { + "_module_label_size_probe.pdf" + } + font_probe_file <- file.path(hc_output_dir, paste0("cluster_heatmap_", mode, font_probe_suffix)) + font_probe_heatmap_parameters <- list( + module_label_preset = "balanced", + module_label_fontsize = cfg$module_label_fontsize, + module_label_pt_size = cfg$module_label_pt_size, + module_box_width_cm = cfg$module_box_width_cm, + gene_count_mode = "text" + ) + hc <- hcocena::hc_plot_cluster_heatmap( + hc, + file_name = basename(font_probe_file), + module_label_preset = font_probe_heatmap_parameters$module_label_preset, + module_label_fontsize = font_probe_heatmap_parameters$module_label_fontsize, + module_label_pt_size = font_probe_heatmap_parameters$module_label_pt_size, + module_box_width_cm = font_probe_heatmap_parameters$module_box_width_cm, + gene_count_mode = font_probe_heatmap_parameters$gene_count_mode, + return_HM = TRUE + ) + plot_variants[[length(plot_variants) + 1L]] <- .hcr_plot_variant( + plot = "Cluster heatmap", + variant = if (isTRUE(cfg$run_module_split)) { + "module-label font-size probe after split" + } else { + "module-label font-size probe" + }, + pdf_file = font_probe_file, + base_dir = output_dir, + stage = if (isTRUE(cfg$run_module_split)) "post module split" else "post clustering", + parameters = font_probe_heatmap_parameters, + defaults = heatmap_plot_defaults + ) + + if (isTRUE(cfg$run_enrichment)) { + enrichment_parameters <- list( + gene_sets = cfg$enrichment_gene_sets, + top = cfg$enrichment_top, + consistent_terms = TRUE, + cluster_columns = FALSE, + store_panel_objects = "never", + heatmap_module_label_fontsize = cfg$module_label_fontsize + ) + hc <- hcocena::hc_functional_enrichment( + hc, + gene_sets = enrichment_parameters$gene_sets, + top = enrichment_parameters$top, + consistent_terms = enrichment_parameters$consistent_terms, + cluster_columns = enrichment_parameters$cluster_columns, + store_panel_objects = enrichment_parameters$store_panel_objects, + heatmap_module_label_fontsize = enrichment_parameters$heatmap_module_label_fontsize + ) + plot_variants[[length(plot_variants) + 1L]] <- .hcr_plot_variant( + plot = "Functional enrichment", + variant = "combined all DBs", + pdf_file = file.path(hc_output_dir, paste0("Enrichment_All_DBs_top_", cfg$enrichment_top, ".pdf")), + base_dir = output_dir, + stage = if (isTRUE(cfg$run_module_split)) "post module split" else "post clustering", + parameters = enrichment_parameters, + defaults = .hcr_enrichment_plot_defaults() + ) + } + + export <- .hcr_export_outputs( + hc, + output_dir, + mode = mode, + font_probe_file = font_probe_file, + plot_variants = plot_variants + ) + visual_report <- .hcr_write_visual_report( + hc = hc, + output_dir = output_dir, + mode = mode, + font_probe_file = font_probe_file, + split_target = split_target, + parameters = c( + list( + mode = mode, + seed = seed, + data_dir = data_dir, + output_dir = output_dir, + reference_dir = reference_dir, + package_source = package_source, + package_version = as.character(utils::packageVersion("hcocena")), + compare_reference = compare_reference, + update_reference = update_reference, + tolerance = tolerance + ), + cfg + ), + plot_variants = plot_variants + ) + export$files[["visual_report"]] <- visual_report + compare <- NULL + if (isTRUE(update_reference)) { + .hcr_update_reference(output_dir, reference_dir) + } + if (isTRUE(compare_reference) && dir.exists(reference_dir) && !isTRUE(update_reference)) { + compare <- .hcr_compare_outputs(output_dir, reference_dir, tolerance = tolerance) + .hcr_write_csv(compare, file.path(output_dir, "reference_comparison.csv")) + } + + manifest <- list( + mode = mode, + created_at = format(Sys.time(), "%Y-%m-%dT%H:%M:%S%z"), + repo_root = repo_root, + data_dir = data_dir, + output_dir = output_dir, + reference_dir = reference_dir, + package_source = package_source, + package_version = as.character(utils::packageVersion("hcocena")), + seed = as.integer(seed), + config = cfg, + split_target = split_target, + files = export$files, + reference_updated = isTRUE(update_reference), + reference_compared = !is.null(compare), + reference_ok = if (is.null(compare)) NA else all(compare$status == "ok") + ) + manifest_path <- .hcr_write_json(manifest, file.path(output_dir, "manifest.json")) + md5_files <- unlist(c(export$files, list(manifest = manifest_path)), use.names = TRUE) + md5 <- data.frame( + file = names(md5_files), + path = as.character(md5_files), + md5 = as.character(tools::md5sum(md5_files)), + stringsAsFactors = FALSE + ) + .hcr_write_csv(md5, file.path(output_dir, "checksums.csv")) + + result <- list( + hc = hc, + manifest = manifest, + output_dir = output_dir, + reference_dir = reference_dir, + comparison = compare, + files = export$files + ) + class(result) <- c("hcocena_realdata_regression", class(result)) + result +} + +.hcr_main <- function() { + args <- .hcr_parse_args() + if (.hcr_bool(args$help, FALSE)) { + cat( + "Usage: Rscript scripts/run_realdata_regression.R [--mode quick|full]\n", + " [--data-dir PATH] [--output-dir PATH] [--reference-dir PATH]\n", + " [--update-reference] [--no-compare-reference] [--load-source true|false]\n", + "\n", + "Environment:\n", + " HCOCENA_REALDATA_DIR Real data directory override.\n", + " HCOCENA_REALDATA_REFERENCE_DIR Parent directory for references.\n", + " HCOCENA_REALDATA_FULL_ENRICHMENT true/false for full-mode enrichment.\n", + sep = "" + ) + return(invisible(NULL)) + } + + mode <- .hcr_scalar(args$mode, Sys.getenv("HCOCENA_REALDATA_MODE", unset = "quick")) + update_reference <- .hcr_bool(args$update_reference, FALSE) + compare_reference <- !.hcr_bool(args$no_compare_reference, FALSE) + load_source <- .hcr_bool(args$load_source, TRUE) + seed <- as.integer(.hcr_scalar(args$seed, Sys.getenv("HCOCENA_REALDATA_SEED", unset = "168575"))) + tolerance <- as.numeric(.hcr_scalar(args$tolerance, "1e-8")) + + res <- run_hcocena_realdata_regression( + mode = mode, + data_dir = args$data_dir, + output_dir = args$output_dir, + reference_dir = args$reference_dir, + update_reference = update_reference, + compare_reference = compare_reference, + load_source = load_source, + seed = seed, + tolerance = tolerance + ) + + cat("Real-data regression completed.\n") + cat("Output: ", res$output_dir, "\n", sep = "") + if (!is.null(res$comparison)) { + cat("Reference comparison:\n") + print(res$comparison, row.names = FALSE) + if (!all(res$comparison$status == "ok")) { + quit(status = 1) + } + } + invisible(res) +} + +if (identical(sys.nframe(), 0L)) { + .hcr_main() +} diff --git a/tests/testthat/test-realdata-regression.R b/tests/testthat/test-realdata-regression.R new file mode 100644 index 0000000..064d3a4 --- /dev/null +++ b/tests/testthat/test-realdata-regression.R @@ -0,0 +1,60 @@ +test_that("real-data regression runner works when explicitly enabled", { + run_realdata <- tolower(Sys.getenv("HCOCENA_RUN_REALDATA", unset = "false")) %in% + c("1", "true", "yes", "y", "on") + skip_if_not(run_realdata, "Set HCOCENA_RUN_REALDATA=true to run local real-data regression tests.") + + repo_root <- normalizePath(test_path("..", ".."), winslash = "/", mustWork = TRUE) + runner <- file.path(repo_root, "scripts", "run_realdata_regression.R") + skip_if_not(file.exists(runner), "Real-data regression runner is unavailable.") + source(runner) + + mode <- Sys.getenv("HCOCENA_REALDATA_MODE", unset = "quick") + data_dir <- Sys.getenv( + "HCOCENA_REALDATA_DIR", + unset = file.path(dirname(repo_root), "data") + ) + skip_if_not(dir.exists(data_dir), paste0("Real-data directory is unavailable: ", data_dir)) + reference_parent <- Sys.getenv("HCOCENA_REALDATA_REFERENCE_DIR", unset = "") + reference_dir <- if (nzchar(reference_parent)) { + file.path(reference_parent, mode) + } else { + file.path(repo_root, "realdata-reference", mode) + } + + out <- run_hcocena_realdata_regression( + mode = mode, + data_dir = data_dir, + output_dir = file.path(repo_root, "realdata-output", paste0("test-", mode)), + reference_dir = reference_dir, + compare_reference = dir.exists(reference_dir), + update_reference = FALSE, + load_source = TRUE + ) + + expect_s3_class(out, "hcocena_realdata_regression") + expect_true(file.exists(file.path(out$output_dir, "manifest.json"))) + expect_true(file.exists(file.path(out$output_dir, "font_probe.csv"))) + expect_true(file.exists(file.path(out$output_dir, "plot_variants.csv"))) + expect_true(file.exists(file.path(out$output_dir, "plot_variant_parameters.csv"))) + expect_true(file.exists(file.path(out$output_dir, paste0("visual_check_report_", mode, ".pdf")))) + + cfg <- .hcr_mode_config(mode) + plot_variants <- utils::read.csv( + file.path(out$output_dir, "plot_variants.csv"), + check.names = FALSE, + stringsAsFactors = FALSE + ) + expect_true(any(plot_variants$plot == "Cluster heatmap" & plot_variants$variant == "baseline before module split")) + if (isTRUE(cfg$run_module_split)) { + expect_true(any(plot_variants$plot == "Cluster heatmap" & grepl("after split", plot_variants$variant))) + } + if (isTRUE(cfg$run_enrichment)) { + expect_true(any(plot_variants$plot == "Functional enrichment" & plot_variants$variant == "combined all DBs")) + } + png_files <- plot_variants$png_file[nzchar(plot_variants$png_file)] + expect_true(all(file.exists(file.path(out$output_dir, png_files)))) + + if (!is.null(out$comparison)) { + expect_true(all(out$comparison$status == "ok")) + } +}) diff --git a/tests/testthat/test-regressions.R b/tests/testthat/test-regressions.R index 29aaab0..bb354b3 100644 --- a/tests/testthat/test-regressions.R +++ b/tests/testthat/test-regressions.R @@ -211,8 +211,8 @@ test_that("regression: heatmap gets subtle smart column gaps only at useful grou test_that("regression: additional heatmap paths reuse layer gap and layer-title logic", { - fun_enrich_path <- test_path("..", "..", "R", "functional_enrichment.R") - knowledge_path <- test_path("..", "..", "R", "plot_enrichment_upstream_network.R") + fun_enrich_path <- test_path("..", "..", "R", "hc_functional_enrichment.R") + knowledge_path <- test_path("..", "..", "R", "hc_plot_enrichment_upstream_network.R") if (!file.exists(fun_enrich_path) || !file.exists(knowledge_path)) { skip("Source files are not available in the installed-package test context.") } @@ -1562,6 +1562,215 @@ test_that("regression: duplicate GFC condition names no longer break cluster sum }) +test_that("split_modules resolves numeric inputs by current heatmap row order", { + resolve_fun <- get(".hc_resolve_modules_for_split", asNamespace("hcocena")) + order_fun <- get(".hc_split_module_order", asNamespace("hcocena")) + + module_label_map <- c(red = "M-red", blue = "M-blue", green = "M-green") + module_order <- order_fun( + cluster_calc = list(heatmap_row_order = c("blue", "green", "red")), + available_colors = c("red", "blue", "green") + ) + + resolved <- resolve_fun( + modules = c(1, 3), + available_colors = c("red", "blue", "green"), + module_label_map = module_label_map, + module_order = module_order + ) + + expect_equal(module_order, c("blue", "green", "red")) + expect_equal(resolved$target_colors, c("blue", "red")) + expect_equal(resolved$resolved_labels, c("M-blue", "M-red")) + expect_equal(resolved$resolution_table$resolved_index, c(1L, 3L)) + expect_equal(resolved$resolution_table$resolved_color, c("blue", "red")) + + bad <- resolve_fun( + modules = 1.5, + available_colors = c("red", "blue", "green"), + module_label_map = module_label_map, + module_order = module_order + ) + expect_equal(bad$resolution_table$status, "not_found") +}) + + +test_that("split_modules maps resolution values per requested module", { + resolution_fun <- get(".hc_split_resolution_by_module", asNamespace("hcocena")) + target_colors <- c("blue", "red", "green", "gold") + resolved_labels <- c("M2", "M3", "M4", "M8") + resolution_table <- data.frame( + input = resolved_labels, + resolved_color = target_colors, + resolved_label = resolved_labels, + status = "ok", + stringsAsFactors = FALSE + ) + + expect_equal( + resolution_fun(0.8, target_colors, resolved_labels, resolution_table), + stats::setNames(rep(0.8, 4), target_colors) + ) + expect_equal( + resolution_fun(c(0.8, 0.8, 0.9, 0.3), target_colors, resolved_labels, resolution_table), + stats::setNames(c(0.8, 0.8, 0.9, 0.3), target_colors) + ) + expect_equal( + resolution_fun( + c(M2 = 0.8, M3 = 0.8, green = 0.9, M8 = 0.3), + target_colors, + resolved_labels, + resolution_table + ), + stats::setNames(c(0.8, 0.8, 0.9, 0.3), target_colors) + ) + + expect_error( + resolution_fun(c(0.8, 0.9), target_colors, resolved_labels, resolution_table), + "one value per resolved module" + ) + expect_error( + resolution_fun(c(M2 = 0.8, M9 = 0.9), target_colors, resolved_labels, resolution_table), + "Unmatched: M9" + ) + expect_error( + resolution_fun(c(M2 = 0.8, M3 = 0.8), target_colors, resolved_labels, resolution_table), + "Missing" + ) +}) + + +test_that("split_modules normalizes label maps and preserves duplicate GFC columns in children", { + normalize_fun <- get(".hc_normalize_module_label_map_for_split", asNamespace("hcocena")) + child_fun <- get(".hc_build_child_cluster_rows", asNamespace("hcocena")) + + expect_equal( + normalize_fun(c(M1 = "red", M2 = "blue"), available_colors = c("red", "blue")), + c(red = "M1", blue = "M2") + ) + + template_row <- data.frame( + clusters = "M1", + gene_no = 4L, + gene_n = "g1,g2,g3,g4", + cluster_included = "yes", + color = "red", + conditions = "", + grp_means = "", + vertexsize = 3, + stringsAsFactors = FALSE + ) + membership <- c(g1 = 1L, g2 = 1L, g3 = 2L, g4 = 2L) + child_rows <- child_fun( + template_row = template_row, + membership = membership, + member_levels = c("1", "2"), + child_colors = c("#111111", "#222222"), + child_labels = c("#111111" = "M1.1", "#222222" = "M1.2"), + gfc_all = data.frame( + T1 = c(1, 3, 10, 12), + T1 = c(5, 7, 20, 22), + Gene = c("g1", "g2", "g3", "g4"), + check.names = FALSE + ) + ) + + expect_identical(child_rows$conditions, c("T1#T1", "T1#T1")) + expect_identical(child_rows$grp_means, c("2,6", "11,21")) +}) + + +test_that("split_modules resolution_test_only without grid does not split", { + hc <- methods::new("HCoCenaExperiment") + + expect_error( + hcocena::hc_split_modules(hc, modules = "M1", resolution_test_only = TRUE), + "requires `resolution_grid`" + ) +}) + + +test_that("split-module labels trigger preserve-existing heatmap numbering", { + has_split_labels <- get(".hc_module_label_map_has_split_labels", asNamespace("hcocena")) + + expect_true(has_split_labels(c(turquoise = "M10.1", red = "M2"))) + expect_false(has_split_labels(c(turquoise = "M10", red = "M2"))) + # Repeated splits append further numeric suffixes and must still be detected. + expect_true(has_split_labels(c(turquoise = "M10.1.2", red = "M2"))) + expect_false(has_split_labels(NULL)) + expect_false(has_split_labels(character())) +}) + + +test_that("cluster heatmap view resolves module order and full GFC scale", { + resolve_modules <- get(".hc_heatmap_view_resolve_modules", asNamespace("hcocena")) + resolve_scale <- get(".hc_heatmap_view_resolve_gfc_scale_limits", asNamespace("hcocena")) + + cluster_calc <- list( + cluster_information = data.frame( + color = c("red", "blue", "green"), + cluster_included = c("yes", "yes", "yes"), + stringsAsFactors = FALSE + ), + heatmap_matrix = matrix( + c(-5, 1, 3, 2, -1, 4), + nrow = 3, + dimnames = list(c("red", "blue", "green"), c("T1", "T2")) + ), + heatmap_row_order = c("blue", "red", "green"), + module_label_map = c(red = "M1", blue = "M2", green = "M3") + ) + + resolved <- resolve_modules( + modules = c("M3", 1, "red"), + cluster_calc = cluster_calc + ) + expect_equal(resolved$target_colors, c("green", "blue", "red")) + expect_equal(resolved$resolved_labels, c("M3", "M2", "M1")) + + hc <- methods::new("HCoCenaExperiment") + hc@integration@cluster <- S4Vectors::SimpleList(cluster_calc) + expect_equal(resolve_scale(hc, cluster_calc), c(-5, 5)) + cluster_calc[["gfc_scale_limits"]] <- c(-3, 3) + expect_equal(resolve_scale(hc, cluster_calc), c(-3, 3)) + expect_true("modules" %in% names(formals(hcocena::hc_plot_cluster_heatmap_view))) + expect_true("write_module_tables" %in% names(formals(get("plot_cluster_heatmap_new", asNamespace("hcocena"))))) +}) + + +test_that("cluster colour palette is unique (merge_clusters keys modules by colour)", { + palette <- hcocena:::get_cluster_colours() + expect_false(any(duplicated(palette))) + expect_true(all(nzchar(palette))) +}) + + +test_that(".hc_with_seed handles NULL / empty / NA seeds without erroring", { + with_seed <- get(".hc_with_seed", asNamespace("hcocena")) + + expect_identical(with_seed(NULL, 1 + 1), 2) + expect_identical(with_seed(integer(0), 1 + 1), 2) + expect_identical(with_seed(NA, 1 + 1), 2) + + # A real seed still makes the draw reproducible. + a <- with_seed(42, stats::runif(1)) + b <- with_seed(42, stats::runif(1)) + expect_identical(a, b) +}) + + +test_that(".hc_longitudinal_group_singletons leaves all-singleton input unchanged", { + group_singletons <- get(".hc_longitudinal_group_singletons", asNamespace("hcocena")) + + ids <- c(a = "1", b = "2", c = "3") + expect_warning( + out <- group_singletons(ids, snn = NULL, group.singletons = TRUE, verbose = FALSE), + NA + ) + expect_identical(out, ids) +}) + + test_that("regression: duplicate heatmap condition labels get layer prefixes and keep axis multiplicity", { display_fun <- get(".hc_gfc_display_col_labels", asNamespace("hcocena")) count_fun <- get(".hc_gfc_display_count_labels", asNamespace("hcocena")) @@ -1635,7 +1844,7 @@ test_that("regression: enrichment-related heatmap legends use standard font sett test_that("regression: enrichment plot body borders use the same thin line style as the main heatmap", { - fun_enrich_path <- test_path("..", "..", "R", "functional_enrichment.R") + fun_enrich_path <- test_path("..", "..", "R", "hc_functional_enrichment.R") if (!file.exists(fun_enrich_path)) { skip("Source files are not available in the installed-package test context.") } @@ -1658,6 +1867,48 @@ test_that("regression: enrichment plot body borders use the same thin line style }) +test_that("regression: enrichment module labels and all-db headers avoid visual overlap", { + fun_enrich_path <- test_path("..", "..", "R", "hc_functional_enrichment.R") + if (!file.exists(fun_enrich_path)) { + skip("Source files are not available in the installed-package test context.") + } + + fun_enrich_src <- paste( + readLines(fun_enrich_path, warn = FALSE), + collapse = "\n" + ) + + expect_true(grepl("module_label_fit <- .hc_module_label_fit_pt", fun_enrich_src, fixed = TRUE)) + expect_true(grepl("module_box_anno <- .hc_module_label_box_annotation", fun_enrich_src, fixed = TRUE)) + expect_true(grepl("label_fontsize_pt = module_label_fit$base_pt_size", fun_enrich_src, fixed = TRUE)) + expect_true(grepl("panel_title_header_offset_mm_all <-", fun_enrich_src, fixed = TRUE)) + expect_true(grepl("panel_title_offset_mm = panel_title_header_offset_mm_all", fun_enrich_src, fixed = TRUE)) +}) + + +test_that("regression: auxiliary heatmap module labels use the shared box-fit guard", { + llm_path <- test_path("..", "..", "R", "hc_plot_module_function_gemini.R") + upstream_path <- test_path("..", "..", "R", "hc_upstream_inference.R") + api_path <- test_path("..", "..", "R", "hc_api.R") + if (!file.exists(llm_path) || !file.exists(upstream_path) || !file.exists(api_path)) { + skip("Source files are not available in the installed-package test context.") + } + + llm_src <- paste(readLines(llm_path, warn = FALSE), collapse = "\n") + upstream_src <- paste(readLines(upstream_path, warn = FALSE), collapse = "\n") + api_src <- paste(readLines(api_path, warn = FALSE), collapse = "\n") + + expect_true(grepl("module_label_fit <- .hc_module_label_fit_pt", llm_src, fixed = TRUE)) + expect_true(grepl("module_box_anno <- .hc_module_label_box_annotation", llm_src, fixed = TRUE)) + expect_true(grepl("label_fontsize_pt = module_label_fit$base_pt_size", llm_src, fixed = TRUE)) + expect_true(grepl("module_label_fit <- .hc_module_label_fit_pt", upstream_src, fixed = TRUE)) + expect_true(grepl("module_box_anno <- .hc_module_label_box_annotation", upstream_src, fixed = TRUE)) + expect_true(grepl("label_fontsize_pt = module_label_fit$base_pt_size", upstream_src, fixed = TRUE)) + expect_true(grepl("is_combined_panel <- panel_key_chr %in% c(\"top_all_dbs\", \"top_all_dbs_mixed\")", api_src, fixed = TRUE)) + expect_true(grepl("max(8, 0.50 * as.numeric(fontsize))", api_src, fixed = TRUE)) +}) + + test_that("regression: lightweight heatmap cache works without ComplexHeatmap object", { cache_info <- get(".hc_heatmap_cache_info", asNamespace("hcocena")) select_col_order <- get(".hc_select_heatmap_col_order", asNamespace("hcocena")) @@ -1791,7 +2042,7 @@ test_that("regression: lightweight heatmap cache works without ComplexHeatmap ob mat_use = matrix(c(-1, 0.5, 1, -0.25), nrow = 2) ) expect_equal(style$module_label_fontsize, 4.2) - expect_equal(style$module_label_pt_size, 0.42) + expect_equal(style$module_label_pt_size, 0.90) expect_equal(style$module_box_width_cm, 0.56) expect_equal(style$cell_size_mm, 4.4) expect_equal(style$gfc_colors, c("#112233", "#f7f7f7", "#cc3311")) @@ -1848,6 +2099,28 @@ test_that("regression: lightweight heatmap cache works without ComplexHeatmap ob save = FALSE ) expect_s3_class(p_heat, "hc_llm_heatmap_plot") + + out_dir <- file.path(tempdir(), paste0("hc_llm_heatmap_export_", Sys.getpid())) + dir.create(out_dir, recursive = TRUE, showWarnings = FALSE) + hc@config@paths <- S4Vectors::DataFrame(dir_output = out_dir) + hc@config@global <- S4Vectors::DataFrame(save_folder = "exports") + p_heat_export <- llm_plot( + hc, + fields = "general_processes", + save = TRUE, + file_stem = "llm_heatmap_test" + ) + heatmap_exports <- attr(p_heat_export, "output_files", exact = TRUE) + expect_true(file.exists(heatmap_exports$pdf)) + expect_true(file.exists(heatmap_exports$png)) + if (requireNamespace("png", quietly = TRUE)) { + read_png <- getExportedValue("png", "readPNG") + png_img <- read_png(heatmap_exports$png) + rgb_img <- png_img[, , seq_len(min(3L, dim(png_img)[[3L]])), drop = FALSE] + expect_gt(stats::sd(as.numeric(rgb_img)), 0.02) + expect_gt(mean(rgb_img), 0.1) + expect_lt(mean(rgb_img), 0.99) + } }) @@ -1882,6 +2155,48 @@ test_that("regression: llm plot export writes pdf and png into configured output expect_true(file.exists(export_files$png)) expect_match(export_files$pdf, "exports") expect_match(export_files$pdf, "llm_test_M1_general_processes\\.pdf$") + + if (requireNamespace("png", quietly = TRUE)) { + read_png <- getExportedValue("png", "readPNG") + png_img <- read_png(export_files$png) + rgb_img <- png_img[, , seq_len(min(3L, dim(png_img)[[3L]])), drop = FALSE] + border_rgb <- c( + as.numeric(rgb_img[1L, , , drop = FALSE]), + as.numeric(rgb_img[dim(rgb_img)[[1L]], , , drop = FALSE]), + as.numeric(rgb_img[, 1L, , drop = FALSE]), + as.numeric(rgb_img[, dim(rgb_img)[[2L]], , drop = FALSE]) + ) + expect_gt(mean(border_rgb), 0.85) + expect_gt(mean(rgb_img), 0.4) + } +}) + + +test_that("regression: llm plot supports RAG comparison fields", { + hc <- methods::new("HCoCenaExperiment") + hc@satellite <- S4Vectors::SimpleList(list( + llm_module_function = list(module_1 = list(status = "ok")), + llm_module_function_summary = data.frame( + module = "M1", + module_color = "steelblue", + general_processes = "Interferon signaling", + contextual_state = "Context-only state", + key_regulators = "STAT1 / IRF7", + with_rag_contextual_state = "RAG-supported contextual state", + stringsAsFactors = FALSE + ) + )) + + plots <- hcocena::hc_plot_module_function_llm( + hc, + with_heatmap = FALSE, + fields = c("contextual_state", "contextual_state_rag"), + save = FALSE + ) + + expect_named(plots, c("contextual_state", "contextual_state_rag")) + expect_true("RAG-supported contextual state" %in% plots$contextual_state_rag$data$term_plot) + expect_true("Context-only state" %in% plots$contextual_state$data$term_plot) }) @@ -2079,6 +2394,77 @@ test_that("regression: split labels and significance suffixes expand module boxe }) +test_that("regression: module label auto-fit uses one safe size for all labels", { + fit_fun <- get(".hc_module_label_fit_pt", asNamespace("hcocena")) + effective_fontsize_fun <- get(".hc_module_label_effective_fontsize", asNamespace("hcocena")) + width_fun <- get(".hc_module_label_text_width_cm", asNamespace("hcocena")) + + fit <- fit_fun( + module_label_pt_size = 0.55, + module_box_width_cm = 0.9, + module_labels_display = c("M1", "M4.12***"), + n_heat_rows = 13, + cell_size_mm = 5.4 + ) + fitted_widths <- width_fun(c("M1", "M4.12***"), fontsize_pt = fit$pt_size) + + expect_true(fit$shrunk) + expect_true(fit$width_limited || fit$height_limited) + expect_equal(fit$pt_size[[1]], fit$pt_size[[2]]) + expect_lt(fit$pt_size[[1]], fit$base_pt_size[[1]]) + expect_true(all(fitted_widths <= fit$available_width_cm + 1e-6)) + expect_lte(max(fit$pt_size), fit$available_height_pt + 1e-6) + + short_fit <- fit_fun( + module_label_pt_size = 0.22, + module_box_width_cm = 0.9, + module_labels_display = c("M1", "M2"), + n_heat_rows = 20, + cell_size_mm = 5.4 + ) + expect_false(short_fit$shrunk) + expect_equal(short_fit$pt_size[[1]], short_fit$pt_size[[2]]) + + fontsize_fit <- fit_fun( + module_label_pt_size = 0.22, + module_box_width_cm = 1.2, + module_labels_display = "M1", + n_heat_rows = 10, + cell_size_mm = 5.4, + module_label_fontsize = 11, + use_fontsize_request = TRUE + ) + expect_equal(fontsize_fit$base_pt_size[[1]], 11) + + fill_fit <- fit_fun( + module_label_pt_size = 0.95, + module_box_width_cm = 0.80, + module_labels_display = c("M1", "M2", "M3"), + n_heat_rows = 5, + cell_size_mm = 7 + ) + fill_widths <- width_fun(c("M1", "M2", "M3"), fontsize_pt = fill_fit$pt_size) + expect_equal(fill_fit$pt_size[[1]], fill_fit$pt_size[[2]]) + expect_gt(fill_fit$pt_size[[1]], 10) + expect_true(all(fill_widths <= fill_fit$available_width_cm + 1e-6)) + expect_equal(effective_fontsize_fun(fill_fit, fallback_fontsize = 5), fill_fit$pt_size[[1]]) + + tiny_box_fit <- fit_fun( + module_label_pt_size = 0.9, + module_box_width_cm = 0.08, + module_labels_display = c("M1", "M123456789"), + n_heat_rows = 6, + cell_size_mm = 5.4, + module_label_fontsize = 80, + use_fontsize_request = TRUE + ) + tiny_box_widths <- width_fun(c("M1", "M123456789"), fontsize_pt = tiny_box_fit$pt_size) + expect_equal(tiny_box_fit$pt_size[[1]], tiny_box_fit$pt_size[[2]]) + expect_true(tiny_box_fit$below_preferred_min) + expect_true(all(tiny_box_widths <= tiny_box_fit$available_width_cm + 1e-6)) +}) + + test_that("regression: current-device heatmap view scales down when device is smaller than export size", { screen_fit_fun <- get(".hc_heatmap_screen_fit_scale", asNamespace("hcocena")) @@ -2242,6 +2628,197 @@ test_that("regression: llm request helpers use ellmer backends", { }) +test_that("regression: llm RAG mode injects retrieved passages and stores citations", { + fun <- get("hc_module_function_llm", asNamespace("hcocena")) + default_rag_url <- get(".hc_llm_default_rag_url", asNamespace("hcocena")) + withr::local_envvar(HCOCENA_RAG_URL = NA_character_) + expect_equal( + default_rag_url(), + "https://limesbcnr-007901.iaas.uni-bonn.de/api/rag/query" + ) + + expect_true(all(c( + "use_rag", + "rag_query", + "rag_url", + "rag_top_k", + "rag_themes", + "rag_min_relevance", + "rag_connect_timeout_sec", + "rag_continue_on_error", + "compare_interpretation_levels" + ) %in% names(formals(fun)))) + + captured_prompts <- character(0) + captured_payload <- NULL + testthat::local_mocked_bindings( + .hc_llm_request_rag = function(query, url, category, top_k, themes, timeout_sec, connect_timeout_sec = NULL) { + captured_payload <<- list( + query = query, + url = url, + category = category, + top_k = top_k, + themes = themes, + timeout_sec = timeout_sec, + connect_timeout_sec = connect_timeout_sec + ) + list( + query = query, + answer = "", + context = list( + list( + chunk = "Prenatal and early postnatal life are key periods of immune development.", + section = "Introduction", + relevance = 0.84, + paper = list( + title = "Neonatal innate immunity", + apa_citation = "Battersby, A. J. (2016). Antimicrobial immunity in early life.", + doi = "10.1000/example" + ) + ), + list( + chunk = "Low relevance passage that should be filtered out.", + section = "Results", + relevance = 0.2, + paper = list( + title = "Unrelated paper", + apa_citation = "Unrelated, A. (2010). Unrelated.", + doi = "" + ) + ) + ), + cited_papers = list( + list( + title = "Neonatal innate immunity", + apa_citation = "Battersby, A. J. (2016). Antimicrobial immunity in early life.", + doi = "10.1000/example" + ) + ) + ) + }, + .hc_llm_request_vllm = function(api_key, + model, + prompt, + system_instruction, + temperature, + timeout_sec, + base_url) { + captured_prompts <<- c(captured_prompts, prompt) + result_text <- if (grepl("Retrieved DoRAG literature context", prompt, fixed = TRUE)) { + '{"general_processes":"rag-supported neonatal innate immune development","contextual_state":"literature-supported interferon-high neonatal state","key_regulators":"STAT1 / IRF7"}' + } else if (grepl("The biological context is: none provided", prompt, fixed = TRUE)) { + '{"general_processes":"gene-only interferon signaling","contextual_state":"interferon-high inflammatory state","key_regulators":"STAT1 / IRF7"}' + } else { + '{"general_processes":"context-aware neonatal immune development","contextual_state":"neonatal interferon-high immune state","key_regulators":"STAT1 / IRF7"}' + } + list(result_text = result_text, raw_response_text = result_text) + }, + .package = "hcocena" + ) + + out <- fun( + genes = c("STAT1", "IRF7", "CXCL10"), + context = "neonatal immune system development", + llm = "vllm", + use_rag = TRUE, + rag_query = "neonatal immune system development", + rag_url = "http://example.test/api/rag/query", + rag_top_k = 2, + rag_themes = "Innate Immunity", + rag_connect_timeout_sec = 45, + rag_min_relevance = 0.5, + compare_interpretation_levels = TRUE, + save_to_hc = FALSE, + verbose = FALSE + ) + + expect_equal(captured_payload$query, "neonatal immune system development") + expect_equal(captured_payload$url, "http://example.test/api/rag/query") + expect_equal(captured_payload$category, "DoRAG") + expect_equal(captured_payload$top_k, 2L) + expect_equal(captured_payload$themes, "Innate Immunity") + expect_equal(captured_payload$connect_timeout_sec, 45) + expect_equal(length(captured_prompts), 3L) + expect_true(grepl("The biological context is: none provided", captured_prompts[[1]], fixed = TRUE)) + expect_false(grepl("neonatal immune system development", captured_prompts[[1]], fixed = TRUE)) + expect_true(grepl("neonatal immune system development", captured_prompts[[2]], fixed = TRUE)) + expect_false(grepl("Retrieved DoRAG literature context", captured_prompts[[2]], fixed = TRUE)) + expect_true(grepl("Retrieved DoRAG literature context", captured_prompts[[3]], fixed = TRUE)) + expect_true(grepl("Prenatal and early postnatal life", captured_prompts[[3]], fixed = TRUE)) + expect_false(grepl("Low relevance passage", captured_prompts[[3]], fixed = TRUE)) + + expect_equal(out$rag_query, "neonatal immune system development") + expect_equal(length(out$rag$context), 1L) + expect_true(grepl("Battersby", out$rag_context_text, fixed = TRUE)) + expect_equal(out$primary_interpretation_level, "with_rag") + expect_named(out$interpretation_levels, c("without_context", "with_context", "with_rag")) + expect_equal(out$response$general_processes, "rag-supported neonatal innate immune development") + expect_equal( + out$interpretation_levels$without_context$response$general_processes, + "gene-only interferon signaling" + ) + expect_equal( + out$interpretation_levels$with_context$response$general_processes, + "context-aware neonatal immune development" + ) + + summary_fun <- get(".hc_llm_summary_from_results", asNamespace("hcocena")) + summary_tbl <- summary_fun(list(out), hc = NULL) + expect_true(summary_tbl$rag_used[[1]]) + expect_equal(summary_tbl$rag_context_count[[1]], 1L) + expect_true(grepl("Battersby", summary_tbl$rag_citations[[1]], fixed = TRUE)) + expect_equal(summary_tbl$primary_interpretation_level[[1]], "with_rag") + expect_equal(summary_tbl$without_context_general_processes[[1]], "gene-only interferon signaling") + expect_equal(summary_tbl$with_context_general_processes[[1]], "context-aware neonatal immune development") + expect_equal(summary_tbl$with_rag_general_processes[[1]], "rag-supported neonatal innate immune development") +}) + + +test_that("regression: llm RAG continue-on-error falls back to context interpretation", { + fun <- get("hc_module_function_llm", asNamespace("hcocena")) + captured_prompts <- character(0) + + testthat::local_mocked_bindings( + .hc_llm_request_rag = function(...) { + stop("DoRAG retrieval request failed: Timeout was reached", call. = FALSE) + }, + .hc_llm_request_vllm = function(api_key, + model, + prompt, + system_instruction, + temperature, + timeout_sec, + base_url) { + captured_prompts <<- c(captured_prompts, prompt) + result_text <- if (grepl("The biological context is: none provided", prompt, fixed = TRUE)) { + '{"general_processes":"gene-only signal","contextual_state":"gene-only state","key_regulators":"STAT1"}' + } else { + '{"general_processes":"context-aware signal","contextual_state":"context-aware state","key_regulators":"STAT1 / IRF7"}' + } + list(result_text = result_text, raw_response_text = result_text) + }, + .package = "hcocena" + ) + + out <- fun( + genes = c("STAT1", "IRF7"), + context = "neonatal immune system development", + llm = "vllm", + use_rag = TRUE, + rag_continue_on_error = TRUE, + compare_interpretation_levels = TRUE, + save_to_hc = FALSE, + verbose = FALSE + ) + + expect_equal(length(captured_prompts), 2L) + expect_equal(out$primary_interpretation_level, "with_context") + expect_named(out$interpretation_levels, c("without_context", "with_context")) + expect_true(grepl("Timeout was reached", out$rag_error_message, fixed = TRUE)) + expect_equal(out$response$general_processes, "context-aware signal") +}) + + test_that("regression: llm summary builder is robust for error-only results", { summary_fun <- get(".hc_llm_summary_from_results", asNamespace("hcocena")) err_res <- list( @@ -2613,3 +3190,43 @@ test_that("regression: htmlwidgets display inline during HTML knitting", { expect_true(any(grepl("hcocena-test-widget html-widget", out, fixed = TRUE))) expect_true(length(knitr::knit_meta()) > 0) }) + +test_that("regression: htmlwidgets dispatch through the print generic for RStudio inline rendering", { + display_src_path <- test_path("..", "..", "R", "hc_display_helpers.R") + if (!file.exists(display_src_path)) { + skip("Source file is not available in the installed-package test context.") + } + display_src <- paste(readLines(display_src_path, warn = FALSE), collapse = "\n") + + # RStudio's notebook renders htmlwidgets inline only when they go through the + # `print` generic (this works even from inside a function). Calling + # print.htmlwidget() directly (e.g. via getFromNamespace) bypasses that hook + # and the widget lands in the Viewer pane instead of inline below the chunk. + expect_true(grepl("print(x)", display_src, fixed = TRUE)) + expect_false(grepl("getFromNamespace(\"print.htmlwidget\"", display_src, fixed = TRUE)) +}) + +test_that("regression: htmlwidgets emit raw HTML only when explicitly opted in", { + skip_if_not_installed("htmlwidgets") + skip_if_not_installed("repr") + + display_fun <- get(".hc_display_object", asNamespace("hcocena")) + widget <- htmlwidgets::createWidget( + name = "hcocena-test-widget", + x = list(value = 1), + package = "htmlwidgets" + ) + + old_options <- options( + rstudio.notebook.executing = TRUE, + viewer = function(...) stop("viewer should not be used", call. = FALSE), + hcocena.htmlwidget_display = "auto", + hcocena.htmlwidget_emit_html = TRUE + ) + on.exit(options(old_options), add = TRUE) + + out <- capture.output(display_fun(widget)) + + expect_true(any(grepl("hcocena-test-widget html-widget", out, fixed = TRUE))) + expect_true(any(grepl("", tolower(out), fixed = TRUE))) +})