Feat@recalc umap (#144)

* Added recalculation of UMAP on manifold tab * Added save and load buttons * Added umap control parameters to import_dataset * bumped version * Added genes_per_anchor parameter
tanaylab · Jul 5, 2023 · 8081ddd · 8081ddd
1 parent 961d781
commit 8081ddd
Show file tree

Hide file tree

Showing 12 changed files with 321 additions and 201 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: MCView
 Title: A Shiny App for Metacell Analysis
-Version: 0.2.12
+Version: 0.2.13
 Authors@R: 
     person(given = "Aviezer",
            family = "Lifshitz",

diff --git a/NEWS.md b/NEWS.md
@@ -1,3 +1,7 @@
+# MCView 0.2.13 
+
+* Added UMAP computation given anchor genes to the 'Manifold' tab.
+
 # MCView 0.2.12
 
 * Added `Sample types` plot to the `Samples` tab.

diff --git a/R/app_server.R b/R/app_server.R
@@ -20,6 +20,8 @@ app_server <- function(input, output, session) {
         globals$screen_height <- input$screen_height
         globals$clipboard <- character(0)
         globals$active_tabs <- config$tabs
+        globals$mc2d <- get_mc_data(dataset(), "mc2d")
+        globals$anchor_genes <- get_mc_data(dataset(), "umap_anchors")
     })
 
     output$menu <- shinydashboard::renderMenu({

diff --git a/R/import.R b/R/import.R
@@ -62,6 +62,10 @@
 #' (in log fraction units - normalized egc) above this threshold
 #' @param minimal_relative_log_fraction When choosing marker genes: take only genes with relative
 #' log fraction (mc_fp) above this this value
+#' @param umap_anchors a vector of gene names to use for UMAP calculation. If NULL, the umap from the anndata object would be used.
+#' @param umap_config a named list with UMAP configuration. See \code{umap::umap} for more details. When NULL, the default configuration would be used, except for: min_dist=0.96, n_neighbors=10, n_epoch=500.
+#' @param min_umap_log_expr minimal log2 expression for genes to use for UMAP calculation.
+#' @param genes_per_anchor number of genes to use for each umap anchor.
 #'
 #' @return invisibly returns an \code{AnnDataR6} object of the read \code{anndata_file}
 #'
@@ -106,6 +110,10 @@ import_dataset <- function(project,
                            copy_atlas = TRUE,
                            minimal_max_log_fraction = -13,
                            minimal_relative_log_fraction = 2,
+                           umap_anchors = NULL,
+                           umap_config = NULL,
+                           min_umap_log_expr = -14,
+                           genes_per_anchor = 30,
                            ...) {
     verbose <- !is.null(getOption("MCView.verbose")) && getOption("MCView.verbose")
     verify_project_dir(project, create = TRUE, atlas = !is.null(atlas_project), ...)
@@ -162,6 +170,8 @@ import_dataset <- function(project,
 
     metacells <- colnames(mc_mat)
 
+    mc_egc <- t(t(mc_mat) / mc_sum)
+
     metadata <- load_metadata(metadata, metadata_fields, metacells, adata)
     if (!is.null(metadata)) {
         serialize_shiny_data(
@@ -183,36 +193,45 @@ import_dataset <- function(project,
     }
 
     cli_alert_info("Processing 2d projection")
-    if (is.null(adata$obsp$obs_outgoing_weights)) {
-        cli_abort_compute_for_mcview("$obsp$obs_outgoing_weights")
+    mc2d_list <- NULL
+    if (!is.null(umap_anchors)) {
+        mc2d_list <- compute_umap(mc_egc, umap_anchors, min_log_expr = min_umap_log_expr, config = umap_config, genes_per_anchor = genes_per_anchor)
+        if (!is.null(mc2d_list)) {
+            serialize_shiny_data(umap_anchors, "umap_anchors", dataset = dataset, cache_dir = cache_dir)
+        }
     }
-    graph <- Matrix::summary(adata$obsp$obs_outgoing_weights) %>%
-        as.data.frame()
-
-    graph <- graph %>% mutate(i = rownames(adata$obs)[i], j = rownames(adata$obs)[j])
 
-    purrr::walk(c("x", "y"), ~ {
-        if (is.null(adata$obs[[.x]])) {
-            cli_abort_compute_for_mcview(glue("$obs${.x}"))
+    if (is.null(mc2d_list)) {
+        if (is.null(adata$obsp$obs_outgoing_weights)) {
+            cli_abort_compute_for_mcview("$obsp$obs_outgoing_weights")
         }
-    })
+        graph <- Matrix::summary(adata$obsp$obs_outgoing_weights) %>%
+            as.data.frame()
 
-    mc2d_list <- list(
-        graph = tibble(mc1 = graph[, 1], mc2 = graph[, 2], weight = graph[, 3]),
-        mc_id = rownames(adata$obs),
-        mc_x = adata$obs %>% select(umap_x = x) %>% tibble::rownames_to_column("mc") %>% tibble::deframe(),
-        mc_y = adata$obs %>% select(umap_y = y) %>% tibble::rownames_to_column("mc") %>% tibble::deframe()
-    )
+        graph <- graph %>% mutate(i = rownames(adata$obs)[i], j = rownames(adata$obs)[j])
+
+        purrr::walk(c("x", "y"), ~ {
+            if (is.null(adata$obs[[.x]])) {
+                cli_abort_compute_for_mcview(glue("$obs${.x}"))
+            }
+        })
+
+        mc2d_list <- list(
+            graph = tibble(mc1 = graph[, 1], mc2 = graph[, 2], weight = graph[, 3]),
+            mc_id = rownames(adata$obs),
+            mc_x = adata$obs %>% select(umap_x = x) %>% tibble::rownames_to_column("mc") %>% tibble::deframe(),
+            mc_y = adata$obs %>% select(umap_y = y) %>% tibble::rownames_to_column("mc") %>% tibble::deframe()
+        )
+    }
     serialize_shiny_data(mc2d_list, "mc2d", dataset = dataset, cache_dir = cache_dir)
 
+
     if (!is.null(metacell_graphs)) {
         cli_alert_info("Processing metacell graphs")
         metacell_graphs <- read_metacell_graphs(metacell_graphs, metacells)
         serialize_shiny_data(metacell_graphs, "metacell_graphs", dataset = dataset, cache_dir = cache_dir)
     }
 
-    mc_egc <- t(t(mc_mat) / mc_sum)
-
     cli_alert_info("Calculating top genes per metacell (marker genes)")
     lateral_field <- "lateral_gene"
 

diff --git a/R/mod_gene_mc.R b/R/mod_gene_mc.R
@@ -67,15 +67,39 @@ mod_gene_mc_server <- function(id, dataset, metacell_types, cell_type_colors, ge
 
             # Projection plots
             output$plot_gene_proj_2d <- render_2d_plotly(input, output, session, dataset, metacell_types, cell_type_colors, gene_modules, globals, source = "proj_mc_plot_gene_tab") %>%
-                bindCache(dataset(), input$color_proj, metacell_types(), cell_type_colors(), input$point_size, input$stroke, input$min_edge_size, input$set_range, input$metacell1, input$metacell2, input$proj_stat, input$expr_range, input$lfp, input$color_proj_gene, input$color_proj_metadata, input$color_proj_gene_module, clipboard_changed(), input$graph_name, input$legend_orientation, input$show_legend_projection, input$scatter_axis_proj, {
-                    if (input$color_proj == "Scatter Axis") {
-                        if (input$scatter_axis_proj == "x") {
-                            c(input$x_axis_var, input$x_axis_type)
-                        } else {
-                            c(input$y_axis_var, input$y_axis_type)
+                bindCache(
+                    dataset(),
+                    input$color_proj,
+                    metacell_types(),
+                    cell_type_colors(),
+                    input$point_size,
+                    input$stroke,
+                    input$min_edge_size,
+                    input$set_range,
+                    input$metacell1,
+                    input$metacell2,
+                    input$proj_stat,
+                    input$expr_range,
+                    input$lfp,
+                    input$color_proj_gene,
+                    input$color_proj_metadata,
+                    input$color_proj_gene_module,
+                    clipboard_changed(),
+                    input$graph_name,
+                    input$legend_orientation,
+                    input$show_legend_projection,
+                    input$scatter_axis_proj,
+                    {
+                        if (input$color_proj == "Scatter Axis") {
+                            if (input$scatter_axis_proj == "x") {
+                                c(input$x_axis_var, input$x_axis_type)
+                            } else {
+                                c(input$y_axis_var, input$y_axis_type)
+                            }
                         }
-                    }
-                })
+                    },
+                    globals$mc2d
+                )
 
             connect_gene_plots(input, output, session, ns, source = "proj_mc_plot_gene_tab")