Merge pull request #1740 from bigomics/devel

New release

Author: ESCRI11

VERSION

@@ -1 +1 @@
-v4.1.4+master260225
+v4.1.4+devel260317

components/app/R/global.R

@@ -25,7 +25,7 @@ if (Sys.info()["sysname"] != "Windows") {
 Sys.setenv("_R_CHECK_LENGTH_1_CONDITION_" = "true")
 
 
-options(shiny.maxRequestSize = 999 * 1024^2) ## max 999Mb upload
+options(shiny.maxRequestSize = 2048 * 1024^2) ## max 2GB (previously 999Mb) upload
 options(shiny.fullstacktrace = TRUE)
 # The following DT global options ensure
 # 1. The header scrolls with the X scroll bar

components/board.biomarker/R/biomarker_server.R

@@ -192,6 +192,13 @@ BiomarkerBoard <- function(id, pgx) {
         )
       }
 
+      shiny::validate(
+        shiny::need(
+          !is.null(res),
+          "Biomarker analysis could not be computed. The current filter settings result in too few features or samples. Please broaden your feature filter or sample selection."
+        )
+      )
+
       is_computed(TRUE)
       return(res)
     })

components/board.clustering/R/clustering_server.R

@@ -141,16 +141,22 @@ ClusteringBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
     # reactive functions ##############
     shiny::observeEvent(
       {
-        list(input$hm_splitby, pgx$X, pgx$samples)
+        list(input$hm_splitby, input$hm_level, pgx$X, pgx$samples)
       },
       {
-        shiny::req(pgx$X, pgx$samples, input$hm_splitby)
+        shiny::req(pgx$X, pgx$samples, input$hm_splitby, input$hm_level)
         if (input$hm_splitby == "none") {
           return()
         }
         if (input$hm_splitby == "gene") {
-          xgenes <- sort(rownames(getFilteredMatrix()$zx))
-          shiny::updateSelectizeInput(session, "hm_splitvar", choices = xgenes, server = TRUE)
+          if (input$hm_level == "geneset") {
+            ## at geneset level, show geneset names instead of genes
+            xgenesets <- sort(rownames(getFilteredMatrix()$zx))
+            shiny::updateSelectizeInput(session, "hm_splitvar", choices = xgenesets, server = TRUE)
+          } else {
+            xgenes <- sort(rownames(getFilteredMatrix()$zx))
+            shiny::updateSelectizeInput(session, "hm_splitvar", choices = xgenes, server = TRUE)
+          }
         }
         if (input$hm_splitby == "phenotype") {
           cvar <- sort(playbase::pgx.getCategoricalPhenotypes(pgx$samples, min.ncat = 2, max.ncat = 999))
@@ -178,6 +184,17 @@ ClusteringBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
       shiny::updateRadioButtons(session, "hm_splitby", selected = "none")
     })
 
+    ## update split radio button label to match current level
+    shiny::observeEvent(input$hm_level, {
+      shiny::req(input$hm_level)
+      level_label <- tspan(input$hm_level, js = FALSE)
+      choices <- c("none", "phenotype", "contrast", "gene")
+      choices_names <- c("none", "phenotype", "contrast", level_label)
+      names(choices) <- choices_names
+      sel <- input$hm_splitby
+      shiny::updateRadioButtons(session, "hm_splitby", choices = choices, selected = sel)
+    })
+
     shiny::observeEvent(pgx, {
       shiny::req(pgx$datatype)
       datatype <- pgx$datatype
@@ -377,8 +394,12 @@ ClusteringBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
       splitby <- input$hm_splitby
       do.split <- splitby != "none"
 
-      if (splitby == "gene" && !splitvar %in% rownames(pgx$X)) {
-        return(NULL)
+      ## validate splitvar exists in the appropriate matrix
+      if (splitby == "gene") {
+        split_rows <- if (input$hm_level == "geneset") rownames(pgx$gsetX) else rownames(pgx$X)
+        if (!splitvar %in% split_rows) {
+          return(NULL)
+        }
       }
       if (splitby == "phenotype" && !splitvar %in% colnames(pgx$samples)) {
         return(NULL)
@@ -402,9 +423,10 @@ ClusteringBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
         shiny::need(is.null(grp) || any(!is.na(grp)), "Selected grouping and filter combination is not valid (no samples left).")
       )
 
-      ## split on gene expression value: hi vs. low
-      if (do.split && splitvar %in% rownames(pgx$X)) {
-        gx <- playbase::imputeMissing(pgx$X, method = "SVD2")
+      ## split on gene/geneset expression value: hi vs. low
+      split_matrix <- if (input$hm_level == "geneset") pgx$gsetX else pgx$X
+      if (do.split && splitvar %in% rownames(split_matrix)) {
+        gx <- playbase::imputeMissing(split_matrix, method = "SVD2")
         gx <- gx[splitvar, colnames(zx)]
 
         ## TODO if this code is revived again, for some datasets
@@ -460,7 +482,7 @@ ClusteringBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
       }
 
       addsplitgene <- function(gg) {
-        if (do.split && splitvar %in% rownames(pgx$X)) {
+        if (do.split && splitvar %in% rownames(split_matrix)) {
           gg <- unique(c(splitvar, gg))
         }
         return(gg)

components/board.expression/R/expression_server.R

@@ -585,7 +585,8 @@ ExpressionBoard <- function(id, pgx, labeltype = shiny::reactive("feature")) {
       if (length(j) == 0) {
         return(NULL)
       } else {
-        gset <- names(which(pgx$GMT[j, ] != 0))
+        gmt_rows <- as.matrix(pgx$GMT[j, , drop = FALSE])
+        gset <- colnames(gmt_rows)[colSums(gmt_rows != 0, na.rm = TRUE) > 0]
         gset <- intersect(gset, rownames(pgx$gsetX))
       }
 

components/board.featuremap/R/featuremap_plot_gset_sig.R

@@ -44,6 +44,10 @@ featuremap_plot_gset_sig_server <- function(id,
       shiny::req(pgx$X)
 
       pos <- pgx$cluster.gsets$pos[["umap2d"]]
+      shiny::validate(shiny::need( ## Custom species has this empty
+        !is.null(pgx$cluster.gsets),
+        "Cluster genesets are not available."
+      ))
       hilight <- NULL
       pheno <- "tissue"
       pheno <- sigvar()

components/board.featuremap/R/featuremap_plot_table_geneset_map.R

@@ -89,6 +89,10 @@ featuremap_plot_table_geneset_map_server <- function(id,
     ns <- session$ns
     sel.row <- 1
     getUMAP <- function() {
+      shiny::validate(shiny::need( ## Custom species has this empty
+        !is.null(pgx$cluster.gsets),
+        "Cluster genesets are not available."
+      ))
       pos <- pgx$cluster.gsets$pos[["umap2d"]]
       colnames(pos) <- c("x", "y")
       pos

components/board.timeseries/R/module_features.R

@@ -172,12 +172,12 @@ TimeSeriesBoard.features_server <- function(id,
       }
 
       if (input$show_statdetails) {
-        sel.p <- grep("^p[.]", colnames(kstats.full))
-        sel.q <- grep("^q[.]", colnames(kstats.full))
+        sel.p <- grep("^p(\\.|$)", colnames(kstats.full))
+        sel.q <- grep("^q(\\.|$)", colnames(kstats.full))
         pq.tables <- kstats.full[, c(sel.p, sel.q), drop = FALSE]
         if (ik %in% names(pgx$gx.meta$meta)) {
-          sel.p <- grep("^p[.]", colnames(ikstats.full))
-          sel.q <- grep("^q[.]", colnames(ikstats.full))
+          sel.p <- grep("^p(\\.|$)", colnames(ikstats.full))
+          sel.q <- grep("^q(\\.|$)", colnames(ikstats.full))
           ik.pq.tables <- ikstats.full[, c(sel.p, sel.q), drop = FALSE]
           colnames(ik.pq.tables) <- paste0(colnames(ik.pq.tables), ".interaction")
           pq.tables <- cbind(pq.tables, ik.pq.tables)

components/board.upload/R/upload_module_computepgx.R

@@ -1165,18 +1165,6 @@ upload_module_computepgx_server <- function(
           sendSuccessMessageToUser = sendSuccessMessageToUser
         )
 
-        path_to_params <- file.path(raw_dir(), "params.RData")
-        saveRDS(params, file = path_to_params)
-
-        # Normalize paths
-        script_path <- normalizePath(file.path(get_opg_root(), "bin", "pgxcreate_op.R"))
-        tmpdir <- normalizePath(raw_dir())
-
-        # Remove global variables
-        try(rm(annot_table))
-        try(rm(custom_geneset))
-
-        # Start the process and store it in the reactive value
         shinyalert::shinyalert(
           title = "Crunching your data!",
           text = stringr::str_squish("Your dataset will be computed in the background.
@@ -1191,6 +1179,17 @@ upload_module_computepgx_server <- function(
           selected = "load-tab"
         )
 
+        path_to_params <- file.path(raw_dir(), "params.RData")
+        saveRDS(params, file = path_to_params)
+
+        # Normalize paths
+        script_path <- normalizePath(file.path(get_opg_root(), "bin", "pgxcreate_op.R"))
+        tmpdir <- normalizePath(raw_dir())
+
+        # Remove global variables
+        try(rm(annot_table))
+        try(rm(custom_geneset))
+
         process_counter(process_counter() + 1)
         dbg("[compute PGX process] : starting processx nr: ", process_counter())
         dbg("[compute PGX process] : process tmpdir = ", tmpdir)

components/board.upload/R/upload_module_normalizationSC.R

@@ -177,7 +177,6 @@ upload_module_normalizationSC_server <- function(id,
         shiny::req(r_samples())
         counts <- r_counts()
         samples <- r_samples()
-
         kk <- intersect(colnames(counts), rownames(samples))
         counts <- counts[, kk, drop = FALSE]
         samples <- samples[kk, , drop = FALSE]
@@ -187,54 +186,22 @@ upload_module_normalizationSC_server <- function(id,
         if (ncol(counts) > cells_trs) {
           dbg("[normalizationSC_server:ds_norm_Counts:] Random sampling of:", cells_trs, "cells.")
           kk <- sample(colnames(counts), cells_trs)
-          counts <- counts[, kk, drop = FALSE]
+          ## Character subsetting on a large dgCMatrix can hit Matrix dispatch issues;
+          ## use integer indices (via match) which always work reliably for sparse matrices.
+          counts <- counts[, match(kk, colnames(counts)), drop = FALSE]
           samples <- samples[kk, , drop = FALSE]
         }
 
         if (input$ref_atlas != "<select>") {
-          dbg("[normalizationSC_server:ds_norm_Counts:] Inferring cell types with Azimuth!")
           dbg("[normalizationSC_server:ds_norm_Counts:] Reference atlas:", input$ref_atlas)
           counts <- as(counts, "dgCMatrix")
-          shiny::withProgress(message = "Inferring cell types with Azimuth...", value = 0.1, {
-            azm <- playbase::pgx.runAzimuth(counts = counts, reference = input$ref_atlas)
-            dbg("[normalizationSC_server:ds_norm_Counts:] Cell types inferred.")
-          })
-
-          if (is.null(azm)) {
-            message_text <- paste(
-              "Please select the Azimuth reference atlas that fits your data.",
-              "NOTE: The uploaded data might NOT be single cell data. Please check your input matrix."
-            )
-            shiny::validate(
-              shiny::need(FALSE, message_text)
-            )
-          }
-
-          if (class(azm) %in% c("matrix", "data.frame")) {
-            kk <- grep("^predicted.*l*2$", colnames(azm))
-            if (any(kk)) {
-              celltype.azimuth <- azm[, kk]
-              if (length(unique(celltype.azimuth)) > 15) {
-                kk <- grep("^predicted.*l*1$", colnames(azm))
-                celltype.azimuth <- azm[, kk]
-              }
-            } else {
-              kk <- grep("^predicted.*subclass*$", colnames(azm))
-              celltype.azimuth <- azm[, kk]
-              if (length(unique(celltype.azimuth)) > 15) {
-                kk <- grep("^predicted.class*$", colnames(azm))
-                celltype.azimuth <- azm[, kk]
-              }
-            }
-            samples <- cbind(samples, celltype.azimuth = celltype.azimuth)
-          } else if (is.vector(azm)) {
-            dbg("[normalizationSC_server:ds_norm_Counts] Azimuth atlas might be incorrect. Please double check.")
-            samples <- cbind(samples, celltype.azimuth = azm)
-          }
 
-          ## Normalization & Dimensional reduction
+          ## logCPM + dim-reduction first: these run fine in the main process and
+          ## their intermediates are freed before the subprocess starts.
           dbg("[normalizationSC_server] Performing logCPM normalization...")
-          nX <- playbase::logCPM(as.matrix(counts), prior = 1, total = 1e4)
+          counts_dense <- as.matrix(counts)
+          nX <- playbase::logCPM(counts_dense, prior = 1, total = 1e4)
+          rm(counts_dense)
           jj <- head(order(-matrixStats::rowSds(nX, na.rm = TRUE)), 250)
           nX1 <- nX[jj, , drop = FALSE]
           nX1 <- nX1 - rowMeans(nX1, na.rm = TRUE)
@@ -247,21 +214,77 @@ upload_module_normalizationSC_server <- function(id,
             pos.list[["umap"]] <- uwot::umap(t(nX1), n_neighbors = max(2, nb))
             pos.list <- lapply(pos.list, function(x) {
               rownames(x) <- colnames(nX1)
-              return(x)
+              x
             })
           })
+          rm(nX, nX1)
+          gc()
           dbg("[normalizationSC_server] PCA, tSNE & UMAP completed.")
 
-          dbg("[normalizationSC_server] Creating & preprocessing Seurat object..")
-          options(Seurat.object.assay.calcn = TRUE)
-          getOption("Seurat.object.assay.calcn")
-          counts <- as(counts, "dgCMatrix")
-          shiny::withProgress(message = "Creating & preprocessing Seurat object...", value = 0.9, {
-            SO <- playbase::pgx.createSeuratObject(counts, samples,
-              batch = NULL, filter = FALSE, preprocess = FALSE
+          ## Run Azimuth + Seurat in a single fresh subprocess.
+          ## Large h5ad uploads leave glibc's heap fragmented after deallocation;
+          ## any subsequent Seurat C++ call (ScaleData, FindNeighbors, RunPCA) then
+          ## crashes. One callr subprocess isolates both Azimuth and the Seurat
+          ## pipeline from the main process heap.
+          dbg("[normalizationSC_server] Running Azimuth + Seurat in subprocess...")
+          shiny::withProgress(message = "Inferring cell types & building Seurat object...", value = 0.7, {
+            result <- callr::r(
+              function(counts, samples, reference, lib) {
+                .libPaths(lib)
+                options(Seurat.object.assay.calcn = TRUE)
+
+                ## Azimuth cell type annotation
+                azm <- tryCatch(
+                  playbase::pgx.runAzimuth(counts = counts, reference = reference),
+                  error = function(e) {
+                    message("[callr] pgx.runAzimuth failed: ", conditionMessage(e))
+                    NULL
+                  }
+                )
+
+                ## Add cell type column to samples metadata
+                if (!is.null(azm)) {
+                  if (is.data.frame(azm) || is.matrix(azm)) {
+                    kk <- grep("^predicted.*l*2$", colnames(azm))
+                    if (any(kk)) {
+                      celltype.azimuth <- azm[, kk]
+                      if (length(unique(celltype.azimuth)) > 15) {
+                        kk <- grep("^predicted.*l*1$", colnames(azm))
+                        celltype.azimuth <- azm[, kk]
+                      }
+                    } else {
+                      kk <- grep("^predicted.*subclass*$", colnames(azm))
+                      celltype.azimuth <- azm[, kk]
+                      if (length(unique(celltype.azimuth)) > 15) {
+                        kk <- grep("^predicted.class*$", colnames(azm))
+                        celltype.azimuth <- azm[, kk]
+                      }
+                    }
+                    samples <- cbind(samples, celltype.azimuth = celltype.azimuth)
+                  } else if (is.vector(azm)) {
+                    samples <- cbind(samples, celltype.azimuth = azm)
+                  }
+                }
+
+                SO <- playbase::pgx.createSeuratObject(counts, samples,
+                  batch = NULL, filter = FALSE, preprocess = FALSE
+                )
+                SO <- playbase::seurat.preprocess(SO, sct = FALSE, tsne = FALSE, umap = FALSE)
+                list(SO = SO, azm_ok = !is.null(azm))
+              },
+              args = list(counts = counts, samples = samples, reference = input$ref_atlas, lib = .libPaths())
             )
-            SO <- playbase::seurat.preprocess(SO, sct = FALSE, tsne = FALSE, umap = FALSE)
           })
+
+          if (!result$azm_ok) {
+            message_text <- paste(
+              "Please select the Azimuth reference atlas that fits your data.",
+              "NOTE: The uploaded data might NOT be single cell data. Please check your input matrix."
+            )
+            shiny::validate(shiny::need(FALSE, message_text))
+          }
+
+          SO <- result$SO
           dbg("[normalizationSC_server] Seurat object created & preprocessed.")
           kk <- setdiff(colnames(samples), colnames(SO@meta.data))
           if (length(kk) > 1) {
@@ -274,7 +297,7 @@ upload_module_normalizationSC_server <- function(id,
             pos.tsne = pos.list[["tsne"]],
             pos.umap = pos.list[["umap"]]
           )
-          rm(counts, nX, nX1, samples, SO)
+          rm(counts, samples, SO)
           return(LL)
         } else {
           shinyalert::shinyalert(
@@ -484,7 +507,7 @@ upload_module_normalizationSC_server <- function(id,
 
       X <- shiny::reactive({
         shiny::req(r_counts())
-        X <- playbase::logCPM(as.matrix(r_counts()), total = 1e4, prior = 1)
+        X <- playbase::logCPM(r_counts(), total = 1e4, prior = 1)
         return(X)
       })