From 30c2e0452acaeb235ad7ec68df2fcc6e2c8e77e6 Mon Sep 17 00:00:00 2001 From: zhanghao-njmu <542370159@qq.com> Date: Wed, 8 Nov 2023 22:05:45 +0800 Subject: [PATCH 1/4] Update the SCExplorer --- R/SCP-analysis.R | 2 +- R/SCP-app.R | 84 ++++++++++++++++++++++-------------------- R/SCP-plot.R | 17 +++++---- R/SCP-workflow.R | 2 + R/utils.R | 4 +- man/FeatureDimPlot.Rd | 2 +- man/FeatureStatPlot.Rd | 3 ++ man/GroupHeatmap.Rd | 2 +- man/RunSCExplorer.Rd | 1 + 9 files changed, 65 insertions(+), 52 deletions(-) diff --git a/R/SCP-analysis.R b/R/SCP-analysis.R index 433a6e0b..28373c27 100644 --- a/R/SCP-analysis.R +++ b/R/SCP-analysis.R @@ -5567,7 +5567,7 @@ RunSCVELO <- function(srt = NULL, assay_X = "RNA", slot_X = "counts", assay_laye stop("One of 'srt', 'adata' must be provided.") } if (is.null(group_by)) { - stop("'roup_by' must be provided.") + stop("'group_by' must be provided.") } if (is.null(linear_reduction) && is.null(nonlinear_reduction)) { stop("'linear_reduction' or 'nonlinear_reduction' must be provided at least one.") diff --git a/R/SCP-app.R b/R/SCP-app.R index 6750d394..b30c405f 100644 --- a/R/SCP-app.R +++ b/R/SCP-app.R @@ -140,11 +140,10 @@ CreateMetaFile <- function(srt, MetaFile, name = NULL, write_tools = FALSE, writ } else { write.attributes <- FALSE } - if (inherits(meta, "numeric")) { + if (is.numeric(meta)) { meta <- as.double(meta) meta_asfeatures <- c(meta_asfeatures, var) - } - if (!is.numeric(meta)) { + } else { if (length(unique(meta)) > ignore_nlevel) { warning("The number of categories in ", var, " is greater than ", ignore_nlevel, ", it will be ignored.", immediate. = TRUE) } else { @@ -249,7 +248,7 @@ PrepareSCExplorer <- function(object, assays = "RNA", slots = c("counts", "data"), ignore_nlevel = 100, write_tools = FALSE, write_misc = FALSE, compression_level = 6, overwrite = FALSE) { - base_dir <- normalizePath(base_dir) + base_dir <- normalizePath(base_dir, mustWork = FALSE) if (!dir.exists(base_dir)) { message("Create SCExplorer base directory: ", base_dir) dir.create(base_dir, recursive = TRUE, showWarnings = FALSE) @@ -585,6 +584,7 @@ RunSCExplorer <- function(base_dir = "SCExplorer", initial_size = 4, initial_ncol = 3, initial_arrange = "Row", + initial_raster = "No", session_workers = 2, plotting_workers = 8, create_script = TRUE, @@ -653,7 +653,6 @@ if (is.null(initial_group)) { if (is.null(initial_feature)) { initial_feature <- meta_features_name[1] } -initial_raster <- ifelse(nrow(data) > 1e5, "Yes", "No") palette_list <- SCP::palette_list theme_list <- list( @@ -892,7 +891,7 @@ ui <- fluidPage( ), fluidRow( column( - width = 6, align = "center", + width = 4, align = "center", radioButtons( inputId = "coExp2", label = "Co-expression", @@ -902,7 +901,17 @@ ui <- fluidPage( ), ), column( - width = 6, align = "center", + width = 4, align = "center", + radioButtons( + inputId = "scale2", + label = "Color scale", + choices = list("feature" = "feature", "all" = "all"), + selected = "feature", + inline = TRUE + ), + ), + column( + width = 4, align = "center", radioButtons( inputId = "raster2", label = "Raster", @@ -1505,70 +1514,70 @@ server <- function(input, output, session) { } } - updateSelectizeInput(session, "features2", choices = c(meta_features_name, all_features), selected = initial_feature, server = TRUE) - updateSelectizeInput(session, "features4", choices = c(meta_features_name, all_features), selected = initial_feature, server = TRUE) - # change dataset ---------------------------------------------------------------- observe({ meta_groups_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset1, "/metadata.stat/asgroups")) reduction_name <- meta_struc[meta_struc$group == paste0("/", input$dataset1, "/reductions"), "name"] default_reduction <- as.character(rhdf5::h5read(MetaFile, name = paste0("/", input$dataset1, "/reductions.stat/Default_reduction"))) all_cells <- rhdf5::h5read(DataFile, name = paste0("/", input$dataset1, "/cells")) - updateSelectInput(session, "reduction1", choices = reduction_name, selected = default_reduction) - updateSelectInput(session, "group1", choices = meta_groups_name, selected = "orig.ident") + updateSelectInput(session, "reduction1", choices = reduction_name, selected = intersect(c(initial_reduction, default_reduction), reduction_name)[1]) + updateSelectInput(session, "group1", choices = meta_groups_name, selected = intersect(c(initial_group, "orig.ident"), meta_groups_name)[1]) updateSelectInput(session, "split1", choices = c("None", meta_groups_name), selected = "None") - updateRadioButtons(session, "raster1", choices = c("Yes", "No"), selected = ifelse(length(all_cells) > 1e5, "Yes", "No")) - }) %>% bindEvent(input$dataset1, ignoreNULL = TRUE, ignoreInit = TRUE) + updateRadioButtons(session, "raster1", choices = c("Yes", "No"), selected = initial_raster) + }) %>% bindEvent(input$dataset1, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ assays <- unique(na.omit(sapply(strsplit(data_group[grep(input$dataset2, data_group)], "/"), function(x) x[3]))) slots <- unique(na.omit(sapply(strsplit(data_group[grep(input$dataset2, data_group)], "/"), function(x) x[4]))) default_assay <- as.character(rhdf5::h5read(DataFile, name = paste0("/", input$dataset2, "/Default_assay"))) default_slot <- ifelse("data" %in% slots, "data", slots[1]) - updateSelectInput(session, "assays2", choices = assays, selected = default_assay) - updateSelectInput(session, "slots2", choices = slots, selected = default_slot) - data <- HDF5Array::TENxMatrix(filepath = DataFile, group = paste0("/", input$dataset2, "/", default_assay, "/", default_slot)) + assay <- intersect(c(initial_assay, default_assay), assays)[1] + slot <- intersect(c(initial_slot, default_slot), slots)[1] + updateSelectInput(session, "assays2", choices = assays, selected = assay) + updateSelectInput(session, "slots2", choices = slots, selected = slot) + data <- HDF5Array::TENxMatrix(filepath = DataFile, group = paste0("/", input$dataset2, "/", assay, "/", slot)) all_features <- colnames(data) all_cells <- rhdf5::h5read(DataFile, name = paste0("/", input$dataset2, "/cells")) meta_features_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset2, "/metadata.stat/asfeatures")) meta_groups_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset2, "/metadata.stat/asgroups")) reduction_name <- meta_struc[meta_struc$group == paste0("/", input$dataset2, "/reductions"), "name"] default_reduction <- as.character(rhdf5::h5read(MetaFile, name = paste0("/", input$dataset2, "/reductions.stat/Default_reduction"))) - updateSelectInput(session, "reduction2", choices = reduction_name, selected = default_reduction) + updateSelectInput(session, "reduction2", choices = reduction_name, selected = intersect(c(initial_reduction, default_reduction), reduction_name)[1]) updateSelectizeInput(session, "features2", - choices = c(meta_features_name, all_features), selected = meta_features_name[1], + choices = c(meta_features_name, all_features), selected = intersect(c(initial_feature, meta_features_name[1]), c(all_features, meta_features_name))[1], options = list(maxOptions = 20, maxItems = 20), server = TRUE ) updateSelectInput(session, "split2", choices = c("None", meta_groups_name), selected = "None") - updateSelectInput(session, "group2", choices = meta_groups_name, selected = "orig.ident") - updateRadioButtons(session, "raster2", choices = c("Yes", "No"), selected = ifelse(length(all_cells) > 1e5, "Yes", "No")) - }) %>% bindEvent(input$dataset2, ignoreNULL = TRUE, ignoreInit = TRUE) + updateRadioButtons(session, "raster2", choices = c("Yes", "No"), selected = initial_raster) + }) %>% bindEvent(input$dataset2, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ meta_groups_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset3, "/metadata.stat/asgroups")) updateSelectInput(session, "stat3", choices = meta_groups_name, selected = "orig.ident") - updateSelectInput(session, "group3", choices = meta_groups_name, selected = "orig.ident") + updateSelectInput(session, "group3", choices = meta_groups_name, selected = intersect(c(initial_group, "orig.ident"), meta_groups_name)[1]) updateSelectInput(session, "split3", choices = c("None", meta_groups_name), selected = "None") - }) %>% bindEvent(input$dataset3, ignoreNULL = TRUE, ignoreInit = TRUE) + }) %>% bindEvent(input$dataset3, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ assays <- unique(na.omit(sapply(strsplit(data_group[grep(input$dataset4, data_group)], "/"), function(x) x[3]))) slots <- unique(na.omit(sapply(strsplit(data_group[grep(input$dataset4, data_group)], "/"), function(x) x[4]))) default_assay <- as.character(rhdf5::h5read(DataFile, name = paste0("/", input$dataset4, "/Default_assay"))) default_slot <- ifelse("data" %in% slots, "data", slots[1]) - updateSelectInput(session, "assays4", choices = assays, selected = default_assay) - updateSelectInput(session, "slots4", choices = slots, selected = default_slot) - data <- HDF5Array::TENxMatrix(filepath = DataFile, group = paste0("/", input$dataset4, "/", default_assay, "/", default_slot)) + assay <- intersect(c(initial_assay, default_assay), assays)[1] + slot <- intersect(c(initial_slot, default_slot), slots)[1] + updateSelectInput(session, "assays4", choices = assays, selected = assay) + updateSelectInput(session, "slots4", choices = slots, selected = slot) + data <- HDF5Array::TENxMatrix(filepath = DataFile, group = paste0("/", input$dataset4, "/", assay, "/", slot)) all_features <- colnames(data) meta_features_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset4, "/metadata.stat/asfeatures")) meta_groups_name <- rhdf5::h5read(MetaFile, name = paste0("/", input$dataset4, "/metadata.stat/asgroups")) updateSelectizeInput(session, "features4", - choices = c(meta_features_name, all_features), selected = meta_features_name[1], + choices = c(meta_features_name, all_features), selected = intersect(c(initial_feature, meta_features_name[1]), c(all_features, meta_features_name))[1], options = list(maxOptions = 20, maxItems = 20), server = TRUE ) updateSelectInput(session, "split4", choices = c("None", meta_groups_name), selected = "None") - updateSelectInput(session, "group4", choices = meta_groups_name, selected = "orig.ident") - }) %>% bindEvent(input$dataset4, ignoreNULL = TRUE, ignoreInit = TRUE) + updateSelectInput(session, "group4", choices = meta_groups_name, selected = intersect(c(initial_group, "orig.ident"), meta_groups_name)[1]) + }) %>% bindEvent(input$dataset4, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ if (input$group3 != "None") { @@ -1747,6 +1756,7 @@ server <- function(input, output, session) { features2 <- input$features2 feature_area2 <- input$feature_area2 coExp2 <- input$coExp2 == "Yes" + scale2 <- input$scale2 raster2 <- input$raster2 == "Yes" palette2 <- input$palette2 theme2 <- input$theme2 @@ -1765,7 +1775,7 @@ server <- function(input, output, session) { features2 <- c(as.character(features2), as.character(feature_area2)) features2 <- unique(features2[features2 %in% c(all_features, meta_features_name)]) if (length(features2) == 0) { - features2 <- meta_features_name[1] + features2 <- intersect(c(initial_feature, meta_features_name[1]), c(all_features, meta_features_name))[1] } promisedData[["p2_dim"]] <- NULL @@ -1788,7 +1798,7 @@ server <- function(input, output, session) { # print(system.time( p2_dim <- SCP::FeatureDimPlot( srt = srt_tmp, features = features2, split.by = split2, reduction = reduction2, slot = "data", raster = raster2, pt.size = pt_size2, - calculate_coexp = coExp2, palette = palette2, theme_use = theme2, + calculate_coexp = coExp2, keep_scale = scale2, palette = palette2, theme_use = theme2, ncol = ncol2, byrow = byrow2, force = TRUE ) # )) @@ -1815,7 +1825,7 @@ server <- function(input, output, session) { bindCache( input$dataset2, input$reduction2, input$split2, input$assays2, input$slots2, input$features2, input$feature_area2, - input$palette2, input$theme2, input$coExp2, input$raster2, + input$palette2, input$theme2, input$coExp2, input$scale2, input$raster2, input$pt_size2, input$size2, input$ncol2, input$arrange2 ) %>% bindEvent(input$submit2, ignoreNULL = FALSE, ignoreInit = FALSE) @@ -2085,7 +2095,7 @@ server <- function(input, output, session) { features4 <- c(as.character(features4), as.character(feature_area4)) features4 <- unique(features4[features4 %in% c(all_features, meta_features_name)]) if (length(features4) == 0) { - features4 <- meta_features_name[1] + features4 <- intersect(c(initial_feature, meta_features_name[1]), c(all_features, meta_features_name))[1] } promisedData[["p4"]] <- NULL @@ -2232,6 +2242,7 @@ server <- function(input, output, session) { "library(promises)", "library(BiocParallel)", "library(ggplot2)", + "library(rlang)", args_code, "plan(multisession, workers = session_workers)", "if (.Platform$OS.type == 'windows') { @@ -2245,11 +2256,6 @@ server <- function(input, output, session) { ) temp <- tempfile("SCExplorer") writeLines(app_code, temp) - wd <- getwd() - on.exit(setwd(wd)) - setwd(base_dir) - source(temp) - setwd(wd) if (isTRUE(create_script)) { app_file <- paste0(base_dir, "/app.R") if (!file.exists(app_file) || isTRUE(overwrite)) { diff --git a/R/SCP-plot.R b/R/SCP-plot.R index b47b7007..ce1d9c5e 100644 --- a/R/SCP-plot.R +++ b/R/SCP-plot.R @@ -1982,7 +1982,7 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli show_stat = ifelse(identical(theme_use, "theme_blank"), FALSE, TRUE), palette = ifelse(isTRUE(compare_features), "Set1", "Spectral"), palcolor = NULL, pt.size = NULL, pt.alpha = 1, bg_cutoff = 0, bg_color = "grey80", - keep_scale = NULL, lower_quantile = 0, upper_quantile = 0.99, lower_cutoff = NULL, upper_cutoff = NULL, + keep_scale = "feature", lower_quantile = 0, upper_quantile = 0.99, lower_cutoff = NULL, upper_cutoff = NULL, add_density = FALSE, density_color = "grey80", density_filled = FALSE, density_filled_palette = "Greys", density_filled_palcolor = NULL, cells.highlight = NULL, cols.highlight = "black", sizes.highlight = 1, alpha.highlight = 1, stroke.highlight = 0.5, calculate_coexp = FALSE, compare_features = FALSE, color_blend_mode = c("blend", "average", "screen", "multiply"), @@ -3343,6 +3343,7 @@ FeatureDimPlot3D <- function(srt, features, reduction = NULL, dims = c(1, 2, 3), #' @param pt.size A numeric value specifying the size of the data points. If NULL, the size is automatically determined. Default is NULL. #' @param pt.alpha A numeric value specifying the transparency of the data points. Default is 1. #' @param jitter.width A numeric value specifying the width of the jitter. Default is 0.5. +#' @param jitter.height A numeric value specifying the height of the jitter. Default is 0.1. #' @param add_trend A logical indicating whether to add a trend line to the plot. Default is FALSE. #' @param trend_color A string specifying the color of the trend line. Default is "black". #' @param trend_linewidth A numeric value specifying the width of the trend line. Default is 1. @@ -3466,7 +3467,7 @@ FeatureStatPlot <- function(srt, stat.by, group.by = NULL, split.by = NULL, bg.b palette = "Paired", palcolor = NULL, alpha = 1, bg_palette = "Paired", bg_palcolor = NULL, bg_alpha = 0.2, add_box = FALSE, box_color = "black", box_width = 0.1, box_ptsize = 2, - add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, + add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, jitter.height = 0.1, add_trend = FALSE, trend_color = "black", trend_linewidth = 1, trend_ptsize = 2, add_stat = c("none", "mean", "median"), stat_color = "black", stat_size = 1, cells.highlight = NULL, cols.highlight = "red", sizes.highlight = 1, alpha.highlight = 1, @@ -3519,7 +3520,7 @@ FeatureStatPlot <- function(srt, stat.by, group.by = NULL, split.by = NULL, bg.b palette = palette, palcolor = palcolor, alpha = alpha, bg_palette = bg_palette, bg_palcolor = bg_palcolor, bg_alpha = bg_alpha, add_box = add_box, box_color = box_color, box_width = box_width, box_ptsize = box_ptsize, - add_point = add_point, pt.color = pt.color, pt.size = pt.size, pt.alpha = pt.alpha, jitter.width = jitter.width, + add_point = add_point, pt.color = pt.color, pt.size = pt.size, pt.alpha = pt.alpha, jitter.width = jitter.width, jitter.height = jitter.height, add_trend = add_trend, trend_color = trend_color, trend_linewidth = trend_linewidth, trend_ptsize = trend_ptsize, add_stat = add_stat, stat_color = stat_color, stat_size = stat_size, cells.highlight = cells.highlight, cols.highlight = cols.highlight, sizes.highlight = sizes.highlight, alpha.highlight = alpha.highlight, @@ -3546,7 +3547,7 @@ FeatureStatPlot <- function(srt, stat.by, group.by = NULL, split.by = NULL, bg.b palette = palette, palcolor = palcolor, alpha = alpha, bg_palette = bg_palette, bg_palcolor = bg_palcolor, bg_alpha = bg_alpha, add_box = add_box, box_color = box_color, box_width = box_width, box_ptsize = box_ptsize, - add_point = add_point, pt.color = pt.color, pt.size = pt.size, pt.alpha = pt.alpha, jitter.width = jitter.width, + add_point = add_point, pt.color = pt.color, pt.size = pt.size, pt.alpha = pt.alpha, jitter.width = jitter.width, jitter.height = jitter.height, add_trend = add_trend, trend_color = trend_color, trend_linewidth = trend_linewidth, trend_ptsize = trend_ptsize, add_stat = add_stat, stat_color = stat_color, stat_size = stat_size, cells.highlight = cells.highlight, cols.highlight = cols.highlight, sizes.highlight = sizes.highlight, alpha.highlight = alpha.highlight, @@ -3662,7 +3663,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp palette = "Paired", palcolor = NULL, alpha = 1, bg_palette = "Paired", bg_palcolor = NULL, bg_alpha = 0.2, add_box = FALSE, box_color = "black", box_width = 0.1, box_ptsize = 2, - add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, + add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, jitter.height = 0.1, add_trend = FALSE, trend_color = "black", trend_linewidth = 1, trend_ptsize = 2, add_stat = c("none", "mean", "median"), stat_color = "black", stat_size = 1, cells.highlight = NULL, cols.highlight = "red", sizes.highlight = 1, alpha.highlight = 1, @@ -4162,7 +4163,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp aes(x = .data[["group.by"]], y = .data[["value"]], linetype = rep(f, nrow(dat)), group = .data[["group.unique"]]), inherit.aes = FALSE, color = pt.color, size = pt.size, alpha = pt.alpha, - position = position_jitterdodge(jitter.width = jitter.width, dodge.width = 0.9, seed = 11), show.legend = FALSE + position = position_jitterdodge(jitter.width = jitter.width, jitter.height = jitter.height, dodge.width = 0.9, seed = 11), show.legend = FALSE )) if (!is.null(cells.highlight)) { cell_df <- subset(p$data, rownames(p$data) %in% cells.highlight) @@ -4170,7 +4171,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp p <- p + geom_point( data = cell_df, aes(x = .data[["group.by"]], y = .data[["value"]], linetype = rep(f, nrow(cell_df)), group = .data[["group.unique"]]), inherit.aes = FALSE, color = cols.highlight, size = sizes.highlight, alpha = alpha.highlight, - position = position_jitterdodge(jitter.width = jitter.width, dodge.width = 0.9, seed = 11), show.legend = FALSE + position = position_jitterdodge(jitter.width = jitter.width, jitter.height = jitter.height, dodge.width = 0.9, seed = 11), show.legend = FALSE ) } } @@ -7854,7 +7855,7 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, db = "GO_BP", TERM2GENE = NULL, TERM2NAME = NULL, minGSSize = 10, maxGSSize = 500, GO_simplify = FALSE, GO_simplify_cutoff = "p.adjust < 0.05", simplify_method = "Wang", simplify_similarityCutoff = 0.7, pvalueCutoff = NULL, padjustCutoff = 0.05, topTerm = 5, show_termid = FALSE, topWord = 20, words_excluded = NULL, - nlabel = 0, features_label = NULL, label_size = 10, label_color = "black", + nlabel = 20, features_label = NULL, label_size = 10, label_color = "black", add_bg = FALSE, bg_alpha = 0.5, add_dot = FALSE, dot_size = unit(8, "mm"), add_reticle = FALSE, reticle_color = "grey", diff --git a/R/SCP-workflow.R b/R/SCP-workflow.R index ea368d83..66165414 100644 --- a/R/SCP-workflow.R +++ b/R/SCP-workflow.R @@ -116,6 +116,8 @@ check_srtList <- function(srtList, batch, assay = NULL, type <- "RNA" } else if (assay_type == "ChromatinAssay") { type <- "Chromatin" + } else { + type <- "Unknown" } } diff --git a/R/utils.R b/R/utils.R index 7d8b6c07..567b01f0 100644 --- a/R/utils.R +++ b/R/utils.R @@ -468,13 +468,13 @@ conda_python <- function(envname = NULL, conda = "auto", all = FALSE) { stop(sprintf(fmt, envname)) } conda_envs <- reticulate::conda_list(conda = conda) - conda_envs <- conda_envs[grep(normalizePath(reticulate:::conda_info(conda = conda)$envs_dirs[1]), x = normalizePath(conda_envs$python), fixed = TRUE), , drop = FALSE] + conda_envs <- conda_envs[grep(normalizePath(reticulate:::conda_info(conda = conda)$envs_dirs[1], mustWork = FALSE), x = normalizePath(conda_envs$python, mustWork = FALSE), fixed = TRUE), , drop = FALSE] env <- conda_envs[conda_envs$name == envname, , drop = FALSE] if (nrow(env) == 0) { stop("conda environment \"", envname, "\" not found") } python <- if (all) env$python else env$python[[1L]] - return(normalizePath(as.character(python))) + return(normalizePath(as.character(python), mustWork = FALSE)) } run_Python <- function(command, envir = .GlobalEnv) { diff --git a/man/FeatureDimPlot.Rd b/man/FeatureDimPlot.Rd index 048a7ab6..6620da4b 100644 --- a/man/FeatureDimPlot.Rd +++ b/man/FeatureDimPlot.Rd @@ -20,7 +20,7 @@ FeatureDimPlot( pt.alpha = 1, bg_cutoff = 0, bg_color = "grey80", - keep_scale = NULL, + keep_scale = "feature", lower_quantile = 0, upper_quantile = 0.99, lower_cutoff = NULL, diff --git a/man/FeatureStatPlot.Rd b/man/FeatureStatPlot.Rd index 66989ceb..fbb0a04d 100644 --- a/man/FeatureStatPlot.Rd +++ b/man/FeatureStatPlot.Rd @@ -33,6 +33,7 @@ FeatureStatPlot( pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, + jitter.height = 0.1, add_trend = FALSE, trend_color = "black", trend_linewidth = 1, @@ -134,6 +135,8 @@ FeatureStatPlot( \item{jitter.width}{A numeric value specifying the width of the jitter. Default is 0.5.} +\item{jitter.height}{A numeric value specifying the height of the jitter. Default is 0.1.} + \item{add_trend}{A logical indicating whether to add a trend line to the plot. Default is FALSE.} \item{trend_color}{A string specifying the color of the trend line. Default is "black".} diff --git a/man/GroupHeatmap.Rd b/man/GroupHeatmap.Rd index 17ec0d2c..d454313f 100644 --- a/man/GroupHeatmap.Rd +++ b/man/GroupHeatmap.Rd @@ -80,7 +80,7 @@ GroupHeatmap( show_termid = FALSE, topWord = 20, words_excluded = NULL, - nlabel = 0, + nlabel = 20, features_label = NULL, label_size = 10, label_color = "black", diff --git a/man/RunSCExplorer.Rd b/man/RunSCExplorer.Rd index ca2401ea..03ca1990 100644 --- a/man/RunSCExplorer.Rd +++ b/man/RunSCExplorer.Rd @@ -22,6 +22,7 @@ RunSCExplorer( initial_size = 4, initial_ncol = 3, initial_arrange = "Row", + initial_raster = "No", session_workers = 2, plotting_workers = 8, create_script = TRUE, From ea08bbb2e0bdbc14ebb5a20200cd2c31f1054a2a Mon Sep 17 00:00:00 2001 From: zhanghao-njmu <542370159@qq.com> Date: Thu, 9 Nov 2023 03:33:06 +0800 Subject: [PATCH 2/4] Fix the "keep_scale" bug in FeatureDimPlot --- R/SCP-plot.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/SCP-plot.R b/R/SCP-plot.R index ce1d9c5e..def37837 100644 --- a/R/SCP-plot.R +++ b/R/SCP-plot.R @@ -2121,6 +2121,8 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli if (!is.numeric(dat_exp) && !inherits(dat_exp, "Matrix")) { stop("'features' must be type of numeric variable.") } + dat_exp[, features][dat_exp[, features] <= bg_cutoff] <- NA + if (length(features) > 50 && !isTRUE(force)) { warning("More than 50 features to be plotted", immediate. = TRUE) answer <- askYesNo("Are you sure to continue?", default = FALSE) @@ -2195,7 +2197,6 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli plist <- lapply(levels(dat_sp[[split.by]]), function(s) { dat <- dat_split[[ifelse(split.by == "All.groups", 1, s)]][, , drop = FALSE] for (f in features) { - dat[, f][dat[, f] <= bg_cutoff] <- NA if (any(is.infinite(dat[, f]))) { dat[, f][which(dat[, f] == max(dat[, f], na.rm = TRUE))] <- max(dat[, f][is.finite(dat[, f])], na.rm = TRUE) dat[, f][which(dat[, f] == min(dat[, f], na.rm = TRUE))] <- min(dat[, f][is.finite(dat[, f])], na.rm = TRUE) @@ -2502,7 +2503,6 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli f <- comb[i, "feature"] s <- comb[i, "split"] dat <- dat_split[[ifelse(split.by == "All.groups", 1, s)]][, c(colnames(dat_use), f), drop = FALSE] - dat[, f][dat[, f] <= bg_cutoff] <- NA if (any(is.infinite(dat[, f]))) { dat[, f][dat[, f] == max(dat[, f], na.rm = TRUE)] <- max(dat[, f][is.finite(dat[, f])], na.rm = TRUE) dat[, f][dat[, f] == min(dat[, f], na.rm = TRUE)] <- min(dat[, f][is.finite(dat[, f])], na.rm = TRUE) From 36969ca24dd71ded4cf496a5226f1f6e70065114 Mon Sep 17 00:00:00 2001 From: zhanghao-njmu <542370159@qq.com> Date: Tue, 14 Nov 2023 15:13:01 +0800 Subject: [PATCH 3/4] Fix some bugs --- DESCRIPTION | 4 +- NAMESPACE | 6 + R/SCP-analysis.R | 396 +++++++++++++++++++++++++++++++++++++--- R/SCP-app.R | 122 +++++++------ R/SCP-cell_annotation.R | 16 +- R/SCP-cellqc.R | 2 +- R/SCP-plot.R | 156 ++++++++-------- R/SCP-workflow.R | 8 +- R/Seurat-function.R | 40 ++-- R/utils.R | 31 ++-- man/CellCorHeatmap.Rd | 10 +- man/CellScoring.Rd | 2 +- man/DynamicHeatmap.Rd | 14 +- man/FeatureHeatmap.Rd | 6 +- man/FeatureStatPlot.Rd | 5 +- man/FetchH5.Rd | 8 +- man/GeneConvert.Rd | 2 +- man/GroupHeatmap.Rd | 29 ++- man/RunSCExplorer.Rd | 12 +- man/as_matrix.Rd | 11 +- src/RcppExports.cpp | 2 +- src/asMatrix.cpp | 23 ++- 22 files changed, 647 insertions(+), 258 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 80d75306..55efbad8 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -21,6 +21,8 @@ Imports: ComplexHeatmap (>= 2.13.0), data.table, dplyr (>= 1.1.0), + future, + future.apply, ggnewscale, ggplot2 (>= 3.4.0), ggrepel, @@ -37,6 +39,7 @@ Imports: mgcv, Matrix, patchwork, + pbapply, png, plotly, proxyC, @@ -68,7 +71,6 @@ Suggests: devtools, DDRTree, e1071, - future, ggpubr, ggridges, ggsignif, diff --git a/NAMESPACE b/NAMESPACE index 477dfd27..066312d5 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -57,6 +57,7 @@ export(FeatureDimPlot3D) export(FeatureHeatmap) export(FeatureStatPlot) export(FetchH5) +export(FindExpressedMarkers) export(GSEAPlot) export(GeneConvert) export(GraphPlot) @@ -206,6 +207,7 @@ importFrom(ComplexHeatmap,width.Legends) importFrom(HDF5Array,TENxMatrix) importFrom(HDF5Array,writeTENxMatrix) importFrom(Matrix,Matrix) +importFrom(Matrix,as.matrix) importFrom(Matrix,colMeans) importFrom(Matrix,colSums) importFrom(Matrix,rowMeans) @@ -228,6 +230,7 @@ importFrom(Seurat,Assays) importFrom(Seurat,AverageExpression) importFrom(Seurat,CaseMatch) importFrom(Seurat,Cells) +importFrom(Seurat,Command) importFrom(Seurat,CreateAssayObject) importFrom(Seurat,CreateDimReducObject) importFrom(Seurat,CreateSeuratObject) @@ -327,6 +330,8 @@ importFrom(dplyr,reframe) importFrom(dplyr,slice_head) importFrom(dplyr,summarise) importFrom(dplyr,summarise_at) +importFrom(future,nbrOfWorkers) +importFrom(future.apply,future_sapply) importFrom(ggforce,geom_mark_ellipse) importFrom(ggforce,geom_mark_hull) importFrom(ggnewscale,new_scale) @@ -511,6 +516,7 @@ importFrom(methods,new) importFrom(methods,slot) importFrom(methods,slotNames) importFrom(patchwork,wrap_plots) +importFrom(pbapply,pbsapply) importFrom(plotly,add_trace) importFrom(plotly,as_widget) importFrom(plotly,layout) diff --git a/R/SCP-analysis.R b/R/SCP-analysis.R index 28373c27..a8b3df7e 100644 --- a/R/SCP-analysis.R +++ b/R/SCP-analysis.R @@ -60,7 +60,7 @@ #' by = list(res$geneID_expand[, "symbol"]), FUN = sum #' ) #' rownames(homologs_counts) <- homologs_counts[, 1] -#' homologs_counts <- as(as.matrix(homologs_counts[, -1]), "dgCMatrix") +#' homologs_counts <- as(as_matrix(homologs_counts[, -1]), "dgCMatrix") #' homologs_counts #' #' @importFrom R.cache loadCache saveCache @@ -470,8 +470,8 @@ searchDatasets <- function(datasets, pattern) { #' @importFrom R.cache loadCache saveCache #' @export CC_GenePrefetch <- function(species = "Homo_sapiens", Ensembl_version = 103, mirror = NULL, max_tries = 5, use_cached_gene = TRUE) { - cc_S_genes <- Seurat::cc.genes.updated.2019$s.genes - cc_G2M_genes <- Seurat::cc.genes.updated.2019$g2m.genes + S <- Seurat::cc.genes.updated.2019$s.genes + G2M <- Seurat::cc.genes.updated.2019$g2m.genes res <- NULL if (species != "Homo_sapiens") { if (isTRUE(use_cached_gene)) { @@ -479,7 +479,7 @@ CC_GenePrefetch <- function(species = "Homo_sapiens", Ensembl_version = 103, mir } if (is.null(res)) { res <- GeneConvert( - geneID = unique(c(cc_S_genes, cc_G2M_genes)), + geneID = unique(c(S, G2M)), geneID_from_IDtype = "symbol", geneID_to_IDtype = "symbol", species_from = "Homo_sapiens", @@ -493,13 +493,13 @@ CC_GenePrefetch <- function(species = "Homo_sapiens", Ensembl_version = 103, mir message("Using cached conversion results for ", species) } genes <- res[["geneID_collapse"]] - cc_S_genes <- unlist(genes[cc_S_genes[cc_S_genes %in% rownames(genes)], "symbol"]) - cc_G2M_genes <- unlist(genes[cc_G2M_genes[cc_G2M_genes %in% rownames(genes)], "symbol"]) + S <- unlist(genes[S[S %in% rownames(genes)], "symbol"]) + G2M <- unlist(genes[G2M[G2M %in% rownames(genes)], "symbol"]) } return(list( res = res, - cc_S_genes = cc_S_genes, - cc_G2M_genes = cc_G2M_genes + S = S, + G2M = G2M )) } @@ -621,14 +621,14 @@ AddModuleScore2 <- function(object, slot = "data", features, pool = NULL, nbin = # return(list(features.use, ctrl.use)) # }, BPPARAM = BPPARAM) # features.scores <- bpaggregate( - # x = as.matrix(assay.data[unlist(lapply(features_collapse, function(x) x[[1]])), ]), + # x = as_matrix(assay.data[unlist(lapply(features_collapse, function(x) x[[1]])), ]), # by = list(unlist(lapply(1:length(features_collapse), function(x) rep(x, length(features_collapse[[x]][[1]]))))), # FUN = mean, # BPPARAM = BPPARAM # ) # features.scores <- features.scores[,-1] # ctrl.scores <- bpaggregate( - # x = as.matrix(assay.data[unlist(lapply(features_collapse, function(x) x[[2]])), ]), + # x = as_matrix(assay.data[unlist(lapply(features_collapse, function(x) x[[2]])), ]), # by = list(unlist(lapply(1:length(features_collapse), function(x) rep(x, length(features_collapse[[x]][[2]]))))), # FUN = mean, # BPPARAM = BPPARAM @@ -672,7 +672,7 @@ AddModuleScore2 <- function(object, slot = "data", features, pool = NULL, nbin = #' ccgenes <- CC_GenePrefetch("Mus_musculus") #' pancreas_sub <- CellScoring( #' srt = pancreas_sub, -#' features = list(S = ccgenes$cc_S_genes, G2M = ccgenes$cc_G2M_genes), +#' features = list(S = ccgenes$S, G2M = ccgenes$G2M), #' method = "Seurat", name = "CC" #' ) #' CellDimPlot(pancreas_sub, "CC_classification") @@ -862,7 +862,7 @@ CellScoring <- function(srt, features = NULL, slot = "data", assay = NULL, split colnames(scores) <- make.names(paste(name, names(features_nm), sep = "_")) } else if (method == "AUCell") { check_R("AUCell") - CellRank <- AUCell::AUCell_buildRankings(as.matrix(GetAssayData(srt_sp, slot = slot, assay = assay)), BPPARAM = BPPARAM, plotStats = FALSE) + CellRank <- AUCell::AUCell_buildRankings(as_matrix(GetAssayData(srt_sp, slot = slot, assay = assay)), BPPARAM = BPPARAM, plotStats = FALSE) cells_AUC <- AUCell::AUCell_calcAUC( geneSets = features, rankings = CellRank, @@ -885,7 +885,7 @@ CellScoring <- function(srt, features = NULL, slot = "data", assay = NULL, split scores_mat <- do.call(rbind, lapply(scores_list, function(x) x[intersect(rownames(x), colnames(srt)), features_used, drop = FALSE])) if (isTRUE(new_assay)) { - srt[[name]] <- CreateAssayObject(counts = t(as.matrix(scores_mat[colnames(srt), , drop = FALSE]))) + srt[[name]] <- CreateAssayObject(counts = t(as_matrix(scores_mat[colnames(srt), , drop = FALSE]))) srt[[name]] <- AddMetaData(object = srt[[name]], metadata = data.frame(termnames = features_nm_used[colnames(scores_mat)])) } else { srt <- AddMetaData(object = srt, metadata = scores_mat) @@ -1077,7 +1077,7 @@ FindConservedMarkers2 <- function(object, grouping.var, ident.1, ident.2 = NULL, object.var <- FetchData(object = object, vars = grouping.var) levels.split <- names(x = sort(x = table(object.var[, 1]))) num.groups <- length(levels.split) - assay <- assay %||% DefaultAssay(srt) + assay <- assay %||% DefaultAssay(object) cells <- list() for (i in 1:num.groups) { @@ -1217,6 +1217,353 @@ FindConservedMarkers2 <- function(object, grouping.var, ident.1, ident.2 = NULL, return(markers.combined) } +#' @examples +#' markers <- FindExpressedMarkers(pancreas_sub, cells.1 = WhichCells(pancreas_sub, expression = Phase == "G2M")) +#' head(markers) +#' FeatureStatPlot(pancreas_sub, rownames(markers)[1], group.by = "Phase", add_point = TRUE) +#' @importFrom Matrix rowSums +#' @importFrom SeuratObject PackageCheck FetchData WhichCells SetIdent Idents +#' @importFrom Seurat FindMarkers FoldChange Command +#' @importFrom pbapply pbsapply +#' @export +FindExpressedMarkers <- function(object, ident.1 = NULL, ident.2 = NULL, cells.1 = NULL, cells.2 = NULL, + features = NULL, assay = NULL, slot = "data", min.expression = 0, + test.use = "wilcox", logfc.threshold = 0.25, base = 2, pseudocount.use = 1, mean.fxn = NULL, fc.name = NULL, + min.pct = 0.1, min.diff.pct = -Inf, max.cells.per.ident = Inf, latent.vars = NULL, only.pos = FALSE, + min.cells.group = 3, min.cells.feature = 3, + norm.method = "LogNormalize", verbose = TRUE, ...) { + ########## FindMarkers.Seurat ########## + + assay <- assay %||% DefaultAssay(object) + data.use <- object[[assay]] + + if (!is.null(cells.1)) { + if (is.null(cells.2)) { + cells.2 <- setdiff(colnames(object), cells.1) + } + } else { + cells.1 <- WhichCells(object = object, idents = ident.1) + cells.2 <- WhichCells(object = object, idents = ident.2) + } + + # fetch latent.vars + if (!is.null(x = latent.vars)) { + latent.vars <- FetchData( + object = object, + vars = latent.vars, + cells = c(cells.1, cells.2) + ) + } + + ########## FindMarkers.Assay ########## + + object <- data.use + + pseudocount.use <- pseudocount.use %||% 1 + data.slot <- ifelse( + test = test.use %in% c("negbinom", "poisson", "DESeq2"), + yes = "counts", + no = slot + ) + data.use <- GetAssayData(object = object, slot = data.slot) + data.use <- data.use[rowSums(data.use) > 0, ] + data.use <- as_matrix(data.use) + data.use[data.use <= min.expression] <- NA + counts <- switch( + EXPR = data.slot, + "scale.data" = GetAssayData(object = object, slot = "counts"), + numeric() + ) + + ########## FoldChange.Assay ########## + features <- features %||% rownames(x = data.use) + slot <- data.slot + + # By default run as if LogNormalize is done + log1pdata.mean.fxn <- function(x) { + return(log(x = rowMeans(x = expm1(x), na.rm = TRUE) + pseudocount.use, base = base)) + } + scaledata.mean.fxn <- function(x) { + return(rowMeans(x, na.rm = TRUE)) + } + counts.mean.fxn <- function(x) { + return(log(x = rowMeans(x = x, na.rm = TRUE) + pseudocount.use, base = base)) + } + if (!is.null(x = norm.method)) { + # For anything apart from log normalization set to rowMeans + if (norm.method != "LogNormalize") { + new.mean.fxn <- counts.mean.fxn + } else { + new.mean.fxn <- switch( + EXPR = slot, + "data" = log1pdata.mean.fxn, + "scale.data" = scaledata.mean.fxn, + "counts" = counts.mean.fxn, + counts.mean.fxn + ) + } + } else { + # If no normalization method is passed use slots to decide mean function + new.mean.fxn <- switch( + EXPR = slot, + "data" = log1pdata.mean.fxn, + "scale.data" = scaledata.mean.fxn, + "counts" = counts.mean.fxn, + log1pdata.mean.fxn + ) + } + mean.fxn <- mean.fxn %||% new.mean.fxn + # Omit the decimal value of e from the column name if base == exp(1) + base.text <- ifelse( + test = base == exp(1), + yes = "", + no = base + ) + fc.name <- fc.name %||% ifelse( + test = slot == "scale.data", + yes = "avg_diff", + no = paste0("avg_log", base.text, "FC") + ) + fc.results <- FoldChange.default( + object = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + features = features, + mean.fxn = mean.fxn, + fc.name = fc.name + ) + + ########## FindMarkers.default ########## + + object <- data.use + slot <- data.slot + + Seurat:::ValidateCellGroups( + object = object, + cells.1 = cells.1, + cells.2 = cells.2, + min.cells.group = min.cells.group + ) + + # reset parameters so no feature filtering is performed + if (test.use %in% "DESeq2") { + features <- rownames(x = object) + min.diff.pct <- -Inf + logfc.threshold <- 0 + } + data <- switch( + EXPR = slot, + "scale.data" = counts, + object + ) + # feature selection (based on percentages) + alpha.min <- pmax(fc.results$pct.1, fc.results$pct.2) + names(x = alpha.min) <- rownames(x = fc.results) + features <- names(x = which(x = alpha.min >= min.pct)) + if (length(x = features) == 0) { + warning("No features pass min.pct threshold; returning empty data.frame") + return(fc.results[features, ]) + } + alpha.diff <- alpha.min - pmin(fc.results$pct.1, fc.results$pct.2) + features <- names( + x = which(x = alpha.min >= min.pct & alpha.diff >= min.diff.pct) + ) + if (length(x = features) == 0) { + warning("No features pass min.diff.pct threshold; returning empty data.frame") + return(fc.results[features, ]) + } + # feature selection (based on logFC) + if (slot != "scale.data") { + total.diff <- fc.results[, 1] # first column is logFC + names(total.diff) <- rownames(fc.results) + features.diff <- if (only.pos) { + names(x = which(x = total.diff >= logfc.threshold)) + } else { + names(x = which(x = abs(x = total.diff) >= logfc.threshold)) + } + features <- intersect(x = features, y = features.diff) + if (length(x = features) == 0) { + warning("No features pass logfc.threshold threshold; returning empty data.frame") + return(fc.results[features, ]) + } + } + # subsample cell groups if they are too large + if (max.cells.per.ident < Inf) { + set.seed(seed = random.seed) + if (length(x = cells.1) > max.cells.per.ident) { + cells.1 <- sample(x = cells.1, size = max.cells.per.ident) + } + if (length(x = cells.2) > max.cells.per.ident) { + cells.2 <- sample(x = cells.2, size = max.cells.per.ident) + } + if (!is.null(x = latent.vars)) { + latent.vars <- latent.vars[c(cells.1, cells.2), , drop = FALSE] + } + } + + de.results <- PerformDE( + object = object, + cells.1 = cells.1, + cells.2 = cells.2, + features = features, + test.use = test.use, + verbose = verbose, + min.cells.feature = min.cells.feature, + latent.vars = latent.vars, + ... + ) + de.results <- cbind(de.results, fc.results[rownames(x = de.results), , drop = FALSE]) + if (only.pos) { + de.results <- de.results[de.results[, 2] > 0, , drop = FALSE] + } + if (test.use %in% "roc") { + de.results <- de.results[order(-de.results$power, -de.results[, 1]), ] + } else { + de.results <- de.results[order(de.results$p_val, -de.results[, 1]), ] + de.results$p_val_adj <- p.adjust( + p = de.results$p_val, + method = "bonferroni", + n = nrow(x = object) + ) + } + + return(de.results) +} + +FoldChange.default <- function(object, cells.1, cells.2, mean.fxn, fc.name, features = NULL, ...) { + features <- features %||% rownames(x = object) + # Calculate percent expressed + thresh.min <- 0 + pct.1 <- round( + x = rowSums(x = object[features, cells.1, drop = FALSE] > thresh.min, na.rm = TRUE) / + length(x = cells.1), + digits = 3 + ) + pct.2 <- round( + x = rowSums(x = object[features, cells.2, drop = FALSE] > thresh.min, na.rm = TRUE) / + length(x = cells.2), + digits = 3 + ) + # Calculate fold change + data.1 <- mean.fxn(object[features, cells.1, drop = FALSE]) + data.2 <- mean.fxn(object[features, cells.2, drop = FALSE]) + fc <- (data.1 - data.2) + fc.results <- as.data.frame(x = cbind(fc, pct.1, pct.2)) + colnames(fc.results) <- c(fc.name, "pct.1", "pct.2") + return(fc.results) +} + + +PerformDE <- function(object, cells.1, cells.2, features, test.use, verbose, min.cells.feature, latent.vars, ...) { + if (!(test.use %in% c("negbinom", "poisson", "MAST", "LR")) && !is.null(x = latent.vars)) { + warning( + "'latent.vars' is only used for the following tests: ", + paste(c("negbinom", "poisson", "MAST", "LR"), collapse = ", "), + call. = FALSE, + immediate. = TRUE + ) + } + data.use <- object[features, c(cells.1, cells.2), drop = FALSE] + data.use <- as_matrix(data.use) + + de.results <- switch( + EXPR = test.use, + "wilcox" = WilcoxDETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + verbose = verbose, + ... + ), + "bimod" = Seurat:::DiffExpTest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + verbose = verbose + ), + "roc" = Seurat:::MarkerTest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + verbose = verbose + ), + "t" = Seurat:::DiffTTest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + verbose = verbose + ), + "negbinom" = Seurat:::GLMDETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + min.cells = min.cells.feature, + latent.vars = latent.vars, + test.use = test.use, + verbose = verbose + ), + "poisson" = Seurat:::GLMDETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + min.cells = min.cells.feature, + latent.vars = latent.vars, + test.use = test.use, + verbose = verbose + ), + "MAST" = Seurat:::MASTDETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + latent.vars = latent.vars, + verbose = verbose, + ... + ), + "DESeq2" = Seurat:::DESeq2DETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + verbose = verbose, + ... + ), + "LR" = Seurat:::LRDETest( + data.use = data.use, + cells.1 = cells.1, + cells.2 = cells.2, + latent.vars = latent.vars, + verbose = verbose + ), + stop("Unknown test: ", test.use) + ) + return(de.results) +} + +#' @importFrom future nbrOfWorkers +#' @importFrom future.apply future_sapply +#' @importFrom pbapply pbsapply +WilcoxDETest <- function(data.use, cells.1, cells.2, verbose = TRUE, ...) { + data.use <- data.use[, c(cells.1, cells.2), drop = FALSE] + j <- seq_len(length.out = length(x = cells.1)) + my.sapply <- ifelse( + test = verbose && nbrOfWorkers() == 1, + yes = pbsapply, + no = future_sapply + ) + check_R("limma") + p_val <- my.sapply( + X = 1:nrow(x = data.use), + FUN = function(x) { + keep <- colnames(data.use)[!is.na(data.use[x, ])] + j <- seq_len(length.out = length(x = intersect(cells.1, keep))) + statistics <- data.use[x, keep] + return(min(2 * min(limma::rankSumTestWithCorrelation(index = j, statistics = statistics)), 1)) + } + ) + + return(data.frame(p_val, row.names = rownames(x = data.use))) +} + + #' Differential gene test #' #' This function utilizes the Seurat package to perform a differential expression (DE) test on gene expression data. @@ -1689,6 +2036,7 @@ RunDEtest <- function(srt, group_by = NULL, group1 = NULL, group2 = NULL, cells1 srt@tools[[paste0("DEtest_", group_by)]][[paste0("ConservedMarkers_", test.use)]] <- data.frame() } } + if (markers_type == "disturbed") { sub_BPPARAM <- SerialParam() bpprogressbar(sub_BPPARAM) <- FALSE @@ -4166,8 +4514,8 @@ orderCells <- function(cds, root_state = NULL, num_paths = NULL, reverse = NULL) num_paths <- 1 } adjusted_S <- t(cds@reducedDimS) - dp <- as.matrix(stats::dist(adjusted_S)) - cellPairwiseDistances(cds) <- as.matrix(stats::dist(adjusted_S)) + dp <- as_matrix(stats::dist(adjusted_S)) + cellPairwiseDistances(cds) <- as_matrix(stats::dist(adjusted_S)) gp <- igraph::graph.adjacency(dp, mode = "undirected", weighted = TRUE) dp_mst <- igraph::minimum.spanning.tree(gp) monocle::minSpanningTree(cds) <- dp_mst @@ -4239,7 +4587,7 @@ project2MST <- function(cds, Projection_Method) { cds <- monocle:::findNearestPointOnMST(cds) closest_vertex <- cds@auxOrderingData[["DDRTree"]]$pr_graph_cell_proj_closest_vertex closest_vertex_names <- colnames(Y)[closest_vertex] - closest_vertex_df <- as.matrix(closest_vertex) + closest_vertex_df <- as_matrix(closest_vertex) row.names(closest_vertex_df) <- row.names(closest_vertex) tip_leaves <- names(which(igraph::degree(dp_mst) == 1)) if (!is.function(Projection_Method)) { @@ -4261,13 +4609,13 @@ project2MST <- function(cds, Projection_Method) { distance <- c(distance, stats::dist(rbind(Z_i, tmp))) } if (!inherits(projection, "matrix")) { - projection <- as.matrix(projection) + projection <- as_matrix(projection) } P[, i] <- projection[which(distance == min(distance))[1], ] } } colnames(P) <- colnames(Z) - dp <- as.matrix(stats::dist(t(P))) + dp <- as_matrix(stats::dist(t(P))) min_dist <- min(dp[dp != 0]) dp <- dp + min_dist diag(dp) <- 0 @@ -4707,9 +5055,9 @@ RunDynamicFeatures <- function(srt, lineages, features = NULL, suffix = lineages t <- srt_sub[[l, drop = TRUE]] t <- t[is.finite(t)] t_ordered <- t[order(t)] - Y_ordered <- as.matrix(Y[features, names(t_ordered), drop = FALSE]) + Y_ordered <- as_matrix(Y[features, names(t_ordered), drop = FALSE]) l_libsize <- Y_libsize[names(t_ordered)] - raw_matrix <- as.matrix(cbind(data.frame(pseudotime = t_ordered), t(Y_ordered))) + raw_matrix <- as_matrix(cbind(data.frame(pseudotime = t_ordered), t(Y_ordered))) # df <- data.frame(x = rowMeans(Y_ordered), y = MatrixGenerics::rowVars(Y_ordered)) # p <- ggplot(df, aes(x = .data[["x"]], y = .data[["y"]])) + @@ -5228,7 +5576,7 @@ adata_to_srt <- function(adata) { next } if (!inherits(obsm, "matrix")) { - obsm <- as.matrix(obsm) + obsm <- as_matrix(obsm) } k <- gsub(pattern = "^X_", replacement = "", x = py_to_r_auto(k)) colnames(obsm) <- paste0(k, "_", seq_len(ncol(obsm))) @@ -5277,7 +5625,7 @@ adata_to_srt <- function(adata) { next } if (!inherits(varm, "matrix")) { - varm <- as.matrix(varm) + varm <- as_matrix(varm) } colnames(varm) <- paste0(py_to_r_auto(k), "_", seq_len(ncol(varm))) rownames(varm) <- py_to_r_auto(adata$var_names$values) @@ -5298,7 +5646,7 @@ adata_to_srt <- function(adata) { next } if (!inherits(varp, "matrix")) { - varp <- as.matrix(varp) + varp <- as_matrix(varp) } colnames(varp) <- py_to_r_auto(adata$var_names$values) rownames(varp) <- py_to_r_auto(adata$var_names$values) diff --git a/R/SCP-app.R b/R/SCP-app.R index b30c405f..d2003679 100644 --- a/R/SCP-app.R +++ b/R/SCP-app.R @@ -323,7 +323,13 @@ PrepareSCExplorer <- function(object, #' data("pancreas_sub") #' pancreas_sub <- Standard_SCP(pancreas_sub) #' PrepareSCExplorer(pancreas_sub, base_dir = "./SCExplorer") -#' srt <- FetchH5(DataFile = "./SCExplorer/Data.hdf5", MetaFile = "./SCExplorer/Meta.hdf5", features = c("Ins1", "Ghrl"), metanames = c("SubCellType", "Phase"), reduction = "UMAP") +#' srt <- FetchH5( +#' DataFile = "./SCExplorer/Data.hdf5", +#' MetaFile = "./SCExplorer/Meta.hdf5", +#' features = c("Ins1", "Ghrl"), +#' metanames = c("SubCellType", "Phase"), +#' reduction = "UMAP" +#' ) #' CellDimPlot(srt, group.by = c("SubCellType", "Phase"), reduction = "UMAP") #' FeatureDimPlot(srt, features = c("Ins1", "Ghrl"), reduction = "UMAP") #' } @@ -506,13 +512,14 @@ CreateSeuratObject2 <- function(counts, project = "SeuratProject", assay = "RNA" #' @param initial_feature A string. The initial feature to be loaded into the app. Default is NULL. #' @param initial_assay A string. The initial assay to be loaded into the app. Default is NULL. #' @param initial_slot A string. The initial slot to be loaded into the app. Default is NULL. -#' @param initial_label A string. Whether to add labels in the initial plot. Default is "No". +#' @param initial_label A string. Whether to add labels in the initial plot. Default is FALSE. #' @param initial_cell_palette A string. The initial color palette for cells. Default is "Paired". #' @param initial_feature_palette A string. The initial color palette for features. Default is "Spectral". #' @param initial_theme A string. The initial theme for plots. Default is "theme_scp". #' @param initial_size A numeric. The initial size of plots. Default is 4. #' @param initial_ncol A numeric. The initial number of columns for arranging plots. Default is 3. -#' @param initial_arrange A string. The initial arrangement of plots. Default is "Row". +#' @param initial_arrange A logical. Whether to use "Row" as the initial arrangement. Default is TRUE. +#' @param initial_raster A logical. Whether to perform rasterization in the initial plot. By default, it is set to automatic, meaning it will be TRUE if the number of cells in the initial datasets exceeds 100,000. #' @param session_workers A numeric. The number of workers for concurrent execution in an asynchronous programming session. Default is 2. #' @param plotting_workers A numeric. The number of threads per worker for parallel plotting. Default is 8. #' @param create_script A logical. Whether to create the SCExplorer app script. Default is TRUE. @@ -577,14 +584,14 @@ RunSCExplorer <- function(base_dir = "SCExplorer", initial_feature = NULL, initial_assay = NULL, initial_slot = NULL, - initial_label = "No", + initial_label = FALSE, initial_cell_palette = "Paired", initial_feature_palette = "Spectral", initial_theme = "theme_scp", initial_size = 4, initial_ncol = 3, - initial_arrange = "Row", - initial_raster = "No", + initial_arrange = NULL, + initial_raster = NULL, session_workers = 2, plotting_workers = 8, create_script = TRUE, @@ -634,6 +641,7 @@ if (!initial_assay %in% assays) { if (!initial_slot %in% slots) { stop("initial_slot is not in the dataset ", initial_slot, " in the DataFile") } +all_cells <- rhdf5::h5read(DataFile, name = paste0("/", initial_dataset, "/cells")) data <- HDF5Array::TENxMatrix(filepath = DataFile, group = paste0("/", initial_dataset, "/", initial_assay, "/", initial_slot)) all_features <- colnames(data) @@ -653,6 +661,12 @@ if (is.null(initial_group)) { if (is.null(initial_feature)) { initial_feature <- meta_features_name[1] } +if (is.null(initial_arrange)) { + initial_arrange <- TRUE +} +if (is.null(initial_raster)) { + initial_raster <- length(all_cells) > 1e5 +} palette_list <- SCP::palette_list theme_list <- list( @@ -714,7 +728,7 @@ ui <- fluidPage( radioButtons( inputId = "label1", label = "Label", - choices = c("Yes", "No"), + choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_label, inline = TRUE ) @@ -724,7 +738,7 @@ ui <- fluidPage( radioButtons( inputId = "raster1", label = "Raster", - choices = c("Yes", "No"), + choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_raster, inline = TRUE ) @@ -737,7 +751,7 @@ ui <- fluidPage( inputId = "pt_size1", label = "Point size", value = 1, - min = 0.5, + min = 0.1, max = 10, step = 0.5, width = "150px" @@ -774,7 +788,7 @@ ui <- fluidPage( radioButtons( inputId = "arrange1", label = "Arrange by", - choices = c("Row", "Column"), + choices = c("Row" = TRUE, "Column" = FALSE), selected = initial_arrange ) ) @@ -895,8 +909,8 @@ ui <- fluidPage( radioButtons( inputId = "coExp2", label = "Co-expression", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ), ), @@ -915,7 +929,7 @@ ui <- fluidPage( radioButtons( inputId = "raster2", label = "Raster", - choices = c("Yes", "No"), + choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_raster, inline = TRUE ) @@ -928,7 +942,7 @@ ui <- fluidPage( inputId = "pt_size2", label = "Point size", value = 1, - min = 0.5, + min = 0.1, max = 10, step = 0.5, width = "150px" @@ -965,7 +979,7 @@ ui <- fluidPage( radioButtons( inputId = "arrange2", label = "Arrange by", - choices = c("Row", "Column"), + choices = c("Row" = TRUE, "Column" = FALSE), selected = initial_arrange ) ) @@ -1080,7 +1094,7 @@ ui <- fluidPage( radioButtons( inputId = "label3", label = "Label", - choices = c("Yes", "No"), + choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_label, inline = TRUE ) @@ -1090,8 +1104,8 @@ ui <- fluidPage( radioButtons( inputId = "flip3", label = "Flip", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1178,7 +1192,7 @@ ui <- fluidPage( radioButtons( inputId = "arrange3", label = "Arrange by", - choices = c("Row", "Column"), + choices = c("Row" = TRUE, "Column" = FALSE), selected = initial_arrange ) ) @@ -1312,8 +1326,8 @@ ui <- fluidPage( radioButtons( inputId = "coExp4", label = "Co-expression", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1322,8 +1336,8 @@ ui <- fluidPage( radioButtons( inputId = "stack4", label = "Stack", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1332,8 +1346,8 @@ ui <- fluidPage( radioButtons( inputId = "flip4", label = "Flip", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ) @@ -1344,8 +1358,8 @@ ui <- fluidPage( radioButtons( inputId = "addbox4", label = "Add box", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1354,8 +1368,8 @@ ui <- fluidPage( radioButtons( inputId = "addpoint4", label = "Add point", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1364,8 +1378,8 @@ ui <- fluidPage( radioButtons( inputId = "addtrend4", label = "Add trend", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ) @@ -1414,8 +1428,8 @@ ui <- fluidPage( radioButtons( inputId = "sameylims4", label = "Same y-axis", - choices = c("Yes", "No"), - selected = "No", + choices = c("Yes" = TRUE, "No" = FALSE), + selected = FALSE, inline = TRUE ) ), @@ -1450,7 +1464,7 @@ ui <- fluidPage( radioButtons( inputId = "arrange4", label = "Arrange by", - choices = c("Row", "Column"), + choices = c("Row" = TRUE, "Column" = FALSE), selected = initial_arrange ) ) @@ -1523,7 +1537,7 @@ server <- function(input, output, session) { updateSelectInput(session, "reduction1", choices = reduction_name, selected = intersect(c(initial_reduction, default_reduction), reduction_name)[1]) updateSelectInput(session, "group1", choices = meta_groups_name, selected = intersect(c(initial_group, "orig.ident"), meta_groups_name)[1]) updateSelectInput(session, "split1", choices = c("None", meta_groups_name), selected = "None") - updateRadioButtons(session, "raster1", choices = c("Yes", "No"), selected = initial_raster) + updateRadioButtons(session, "raster1", choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_raster) }) %>% bindEvent(input$dataset1, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ @@ -1548,7 +1562,7 @@ server <- function(input, output, session) { options = list(maxOptions = 20, maxItems = 20), server = TRUE ) updateSelectInput(session, "split2", choices = c("None", meta_groups_name), selected = "None") - updateRadioButtons(session, "raster2", choices = c("Yes", "No"), selected = initial_raster) + updateRadioButtons(session, "raster2", choices = c("Yes" = TRUE, "No" = FALSE), selected = initial_raster) }) %>% bindEvent(input$dataset2, ignoreNULL = TRUE, ignoreInit = FALSE) observe({ @@ -1615,14 +1629,14 @@ server <- function(input, output, session) { } else { split1 <- input$split1 } - label1 <- input$label1 == "Yes" - raster1 <- input$raster1 == "Yes" + label1 <- input$label1 + raster1 <- input$raster1 palette1 <- input$palette1 theme1 <- input$theme1 size1 <- input$size1 pt_size1 <- input$pt_size1 ncol1 <- input$ncol1 - byrow1 <- input$arrange1 == "Row" + byrow1 <- input$arrange1 # lapply(grep("1$",names(input),value = TRUE), function(x)print(paste0(x,":",input[[x]]))) @@ -1755,15 +1769,15 @@ server <- function(input, output, session) { slots2 <- input$slots2 features2 <- input$features2 feature_area2 <- input$feature_area2 - coExp2 <- input$coExp2 == "Yes" + coExp2 <- input$coExp2 scale2 <- input$scale2 - raster2 <- input$raster2 == "Yes" + raster2 <- input$raster2 palette2 <- input$palette2 theme2 <- input$theme2 size2 <- input$size2 pt_size2 <- input$pt_size2 ncol2 <- input$ncol2 - byrow2 <- input$arrange2 == "Row" + byrow2 <- input$arrange2 # lapply(grep("2$",names(input),value = TRUE), function(x)print(paste0(x,":",input[[x]]))) @@ -1812,7 +1826,7 @@ server <- function(input, output, session) { if (isTRUE(plot3d)) { p2_3d <- SCP::FeatureDimPlot3D( srt = srt_tmp, features = features2, reduction = reduction2, pt.size = pt_size2 * 2, - calculate_coexp = ifelse(coExp2 == "Yes", TRUE, FALSE), force = TRUE + calculate_coexp = coExp2, force = TRUE ) } else { p2_3d <- NULL @@ -1918,14 +1932,14 @@ server <- function(input, output, session) { } else { split3 <- input$split3 } - label3 <- input$label3 == "Yes" - flip3 <- input$flip3 == "Yes" + label3 <- input$label3 + flip3 <- input$flip3 palette3 <- input$palette3 theme3 <- input$theme3 labelsize3 <- input$labelsize3 size3 <- input$size3 ncol3 <- input$ncol3 - byrow3 <- input$arrange3 == "Row" + byrow3 <- input$arrange3 if (input$aspectratio3 == "auto") { aspect.ratio <- NULL } else { @@ -2067,18 +2081,18 @@ server <- function(input, output, session) { feature_area4 <- input$feature_area4 plotby4 <- input$plotby4 fillby4 <- input$fillby4 - coExp4 <- input$coExp4 == "Yes" - stack4 <- input$stack4 == "Yes" - flip4 <- input$flip4 == "Yes" - addbox4 <- input$addbox4 == "Yes" - addpoint4 <- input$addpoint4 == "Yes" - addtrend4 <- input$addtrend4 == "Yes" + coExp4 <- input$coExp4 + stack4 <- input$stack4 + flip4 <- input$flip4 + addbox4 <- input$addbox4 + addpoint4 <- input$addpoint4 + addtrend4 <- input$addtrend4 palette4 <- input$palette4 theme4 <- input$theme4 - sameylims4 <- input$sameylims4 == "Yes" + sameylims4 <- input$sameylims4 size4 <- input$size4 ncol4 <- input$ncol4 - byrow4 <- input$arrange4 == "Row" + byrow4 <- input$arrange4 if (input$aspectratio4 == "auto") { aspect.ratio <- NULL } else { diff --git a/R/SCP-cell_annotation.R b/R/SCP-cell_annotation.R index bfef9bff..4ac82383 100644 --- a/R/SCP-cell_annotation.R +++ b/R/SCP-cell_annotation.R @@ -276,7 +276,7 @@ RunKNNPredict <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, } if (!inherits(ref, "matrix")) { - ref <- as.matrix(ref) + ref <- as_matrix(ref) } k <- min(c(k, nrow(ref))) @@ -321,7 +321,7 @@ RunKNNPredict <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, # rownames(tst.expr) <- rownames(ref) # colnames(tst.expr) <- colnames(query) # features_common <- intersect(rownames(query), rownames(ref.expr)) - # tst.expr[features_common, ] <- as.matrix(query[features_common, ]) + # tst.expr[features_common, ] <- as_matrix(query[features_common, ]) # # library(RcppArmadillo) # library(RcppXPtrUtils) @@ -335,7 +335,7 @@ RunKNNPredict <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, # }", # depends = c("RcppArmadillo") # ) - # cors <- as.matrix(parDist(t(query), method = "custom", func = CosineCPP)) + # cors <- as_matrix(parDist(t(query), method = "custom", func = CosineCPP)) # cors <- cors[colnames(ref.expr), colnames(tst.expr)] # @@ -407,11 +407,11 @@ RunKNNPredict <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, ) } if (k == 1) { - match_k_cell <- as.matrix(apply(d, 2, function(x) names(x)[order(x, decreasing = FALSE)[1]])) - match_k_distance <- as.matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1]])) + match_k_cell <- as_matrix(apply(d, 2, function(x) names(x)[order(x, decreasing = FALSE)[1]])) + match_k_distance <- as_matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1]])) } else { - match_k_cell <- t(as.matrix(apply(d, 2, function(x) names(x)[order(x, decreasing = FALSE)[1:k]]))) - match_k_distance <- t(as.matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1:k]]))) + match_k_cell <- t(as_matrix(apply(d, 2, function(x) names(x)[order(x, decreasing = FALSE)[1:k]]))) + match_k_distance <- t(as_matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1:k]]))) } } @@ -438,7 +438,7 @@ RunKNNPredict <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, x <- x / sum(x) return(x) })) - match_prob <- as.matrix(match_prob) + match_prob <- as_matrix(match_prob) rownames(match_prob) <- names(match_freq) match_best <- apply(match_prob, 1, function(x) names(x)[order(x, decreasing = TRUE)][1]) } diff --git a/R/SCP-cellqc.R b/R/SCP-cellqc.R index 00036637..c716d3b8 100644 --- a/R/SCP-cellqc.R +++ b/R/SCP-cellqc.R @@ -96,7 +96,7 @@ db_Scrublet <- function(srt, assay = "RNA", db_rate = ncol(srt) / 1000 * 0.01, . } check_Python("scrublet") scr <- import("scrublet") - raw_counts <- t(as.matrix(GetAssayData(object = srt, assay = assay, slot = "counts"))) + raw_counts <- t(as_matrix(GetAssayData(object = srt, assay = assay, slot = "counts"))) scrub <- scr$Scrublet(raw_counts, expected_doublet_rate = db_rate, ...) res <- scrub$scrub_doublets() doublet_scores <- res[[1]] diff --git a/R/SCP-plot.R b/R/SCP-plot.R index def37837..602d2a61 100644 --- a/R/SCP-plot.R +++ b/R/SCP-plot.R @@ -1548,7 +1548,7 @@ CellDimPlot <- function(srt, group.by, reduction = NULL, dims = c(1, 2), split.b density <- NULL } if (!is.null(graph)) { - net_mat <- as.matrix(x = graph)[rownames(dat), rownames(dat)] + net_mat <- as_matrix(graph)[rownames(dat), rownames(dat)] net_mat[net_mat == 0] <- NA net_mat[upper.tri(net_mat)] <- NA net_df <- reshape2::melt(net_mat, na.rm = TRUE, stringsAsFactors = FALSE) @@ -2098,24 +2098,24 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli if (length(features_gene) > 0) { if (all(rownames(srt@assays[[assay]]) %in% features_gene)) { - dat_gene <- t(slot(srt@assays[[assay]], slot)) + dat_gene <- t(as_matrix(slot(srt@assays[[assay]], slot))) } else { - dat_gene <- t(slot(srt@assays[[assay]], slot)[features_gene, , drop = FALSE]) + dat_gene <- t(as_matrix(slot(srt@assays[[assay]], slot)[features_gene, , drop = FALSE])) } } else { dat_gene <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_meta) > 0) { - dat_meta <- as.matrix(srt@meta.data[, features_meta, drop = FALSE]) + dat_meta <- as_matrix(srt@meta.data[, features_meta, drop = FALSE]) } else { dat_meta <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_embedding) > 0) { - dat_embedding <- as.matrix(FetchData(srt, vars = features_embedding)) + dat_embedding <- as_matrix(FetchData(srt, vars = features_embedding)) } else { dat_embedding <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } - dat_exp <- do.call(cbind, list(dat_gene, dat_meta, dat_embedding)) + dat_exp <- as_matrix(do.call(cbind, list(dat_gene, dat_meta, dat_embedding))) features <- unique(features[features %in% c(features_gene, features_meta, features_embedding)]) if (!is.numeric(dat_exp) && !inherits(dat_exp, "Matrix")) { @@ -2263,7 +2263,7 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli cells.highlight_use <- rownames(dat)[dat[["color_blend"]] != bg_color] } if (!is.null(graph)) { - net_mat <- as.matrix(x = graph)[rownames(dat), rownames(dat)] + net_mat <- as_matrix(graph)[rownames(dat), rownames(dat)] net_mat[net_mat == 0] <- NA net_mat[upper.tri(net_mat)] <- NA net_df <- melt(net_mat, na.rm = TRUE, stringsAsFactors = FALSE) @@ -2534,7 +2534,7 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli colors_value <- seq(lower_cutoff %||% quantile(dat_exp[is.finite(dat_exp[, f]), f], lower_quantile, na.rm = TRUE), upper_cutoff %||% quantile(dat_exp[is.finite(dat_exp[, f]), f], upper_quantile, na.rm = TRUE) + 0.001, length.out = 100) } if (keep_scale == "all") { - all_values <- as.matrix(dat_exp[, features]) + all_values <- as_matrix(dat_exp[, features]) colors_value <- seq(lower_cutoff %||% quantile(all_values[is.finite(all_values)], lower_quantile, na.rm = TRUE), upper_cutoff %||% quantile(all_values, upper_quantile, na.rm = TRUE) + 0.001, length.out = 100) } } @@ -2542,7 +2542,7 @@ FeatureDimPlot <- function(srt, features, reduction = NULL, dims = c(1, 2), spli dat[which(dat[, "value"] > max(colors_value, na.rm = TRUE)), "value"] <- max(colors_value, na.rm = TRUE) dat[which(dat[, "value"] < min(colors_value, na.rm = TRUE)), "value"] <- min(colors_value, na.rm = TRUE) if (!is.null(graph)) { - net_mat <- as.matrix(x = graph)[rownames(dat), rownames(dat)] + net_mat <- as_matrix(graph)[rownames(dat), rownames(dat)] net_mat[net_mat == 0] <- NA net_mat[upper.tri(net_mat)] <- NA net_df <- melt(net_mat, na.rm = TRUE, stringsAsFactors = FALSE) @@ -3002,6 +3002,7 @@ CellDimPlot3D <- function(srt, group.by, reduction = NULL, dims = c(1, 2, 3), ax #' FeatureDimPlot3D(pancreas_sub, features = c("StandardPC_1", "StandardPC_2"), reduction = "StandardpcaUMAP3D") #' #' @importFrom Seurat Reductions Embeddings Key FetchData +#' @importFrom Matrix t #' @importFrom methods slot #' @importFrom utils askYesNo #' @importFrom plotly plot_ly add_trace layout as_widget @@ -3111,24 +3112,24 @@ FeatureDimPlot3D <- function(srt, features, reduction = NULL, dims = c(1, 2, 3), if (length(features_gene) > 0) { if (all(rownames(srt@assays[[assay]]) %in% features_gene)) { - dat_gene <- t(slot(srt@assays[[assay]], slot)) + dat_gene <- t(as_matrix(slot(srt@assays[[assay]], slot))) } else { - dat_gene <- t(slot(srt@assays[[assay]], slot)[features_gene, , drop = FALSE]) + dat_gene <- t(as_matrix(slot(srt@assays[[assay]], slot)[features_gene, , drop = FALSE])) } } else { dat_gene <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_meta) > 0) { - dat_meta <- as.matrix(srt@meta.data[, features_meta, drop = FALSE]) + dat_meta <- as_matrix(srt@meta.data[, features_meta, drop = FALSE]) } else { dat_meta <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_embedding) > 0) { - dat_embedding <- as.matrix(FetchData(srt, vars = features_embedding)) + dat_embedding <- as_matrix(FetchData(srt, vars = features_embedding)) } else { dat_embedding <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } - dat_exp <- do.call(cbind, list(dat_gene, dat_meta, dat_embedding)) + dat_exp <- as_matrix(do.call(cbind, list(dat_gene, dat_meta, dat_embedding))) features <- unique(features[features %in% c(features_gene, features_meta, features_embedding)]) if (!is.numeric(dat_exp) && !inherits(dat_exp, "Matrix")) { @@ -3429,6 +3430,7 @@ FeatureDimPlot3D <- function(srt, features, reduction = NULL, dims = c(1, 2, 3), #' group.by = "SubCellType", bg.by = "CellType", stack = TRUE, flip = TRUE #' ) %>% panel_fix_overall(width = 8, height = 5) # As the plot is created by combining, we can adjust the overall height and width directly. #' +#' FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "group") #' FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "feature") #' FeatureStatPlot(pancreas_sub, #' stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "feature", @@ -3449,7 +3451,7 @@ FeatureDimPlot3D <- function(srt, features, reduction = NULL, dims = c(1, 2, 3), #' #' library(Matrix) #' data <- pancreas_sub@assays$RNA@data -#' pancreas_sub@assays$RNA@scale.data <- as.matrix(data / rowMeans(data)) +#' pancreas_sub@assays$RNA@scale.data <- as_matrix(data / rowMeans(data)) #' FeatureStatPlot(pancreas_sub, #' stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", #' slot = "scale.data", ylab = "FoldChange", same.y.lims = TRUE, y.max = 4 @@ -3467,7 +3469,7 @@ FeatureStatPlot <- function(srt, stat.by, group.by = NULL, split.by = NULL, bg.b palette = "Paired", palcolor = NULL, alpha = 1, bg_palette = "Paired", bg_palcolor = NULL, bg_alpha = 0.2, add_box = FALSE, box_color = "black", box_width = 0.1, box_ptsize = 2, - add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, jitter.height = 0.1, + add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.4, jitter.height = 0.1, add_trend = FALSE, trend_color = "black", trend_linewidth = 1, trend_ptsize = 2, add_stat = c("none", "mean", "median"), stat_color = "black", stat_size = 1, cells.highlight = NULL, cols.highlight = "red", sizes.highlight = 1, alpha.highlight = 1, @@ -3663,7 +3665,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp palette = "Paired", palcolor = NULL, alpha = 1, bg_palette = "Paired", bg_palcolor = NULL, bg_alpha = 0.2, add_box = FALSE, box_color = "black", box_width = 0.1, box_ptsize = 2, - add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.5, jitter.height = 0.1, + add_point = FALSE, pt.color = "grey30", pt.size = NULL, pt.alpha = 1, jitter.width = 0.4, jitter.height = 0.1, add_trend = FALSE, trend_color = "black", trend_linewidth = 1, trend_ptsize = 2, add_stat = c("none", "mean", "median"), stat_color = "black", stat_size = 1, cells.highlight = NULL, cols.highlight = "red", sizes.highlight = 1, alpha.highlight = 1, @@ -3811,7 +3813,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp dat_gene <- matrix(nrow = length(allcells), ncol = 0) } if (length(features_meta) > 0) { - dat_meta <- as.matrix(meta.data[, features_meta, drop = FALSE]) + dat_meta <- as_matrix(meta.data[, features_meta, drop = FALSE]) } else { dat_meta <- matrix(nrow = length(allcells), ncol = 0) } @@ -3846,7 +3848,7 @@ ExpressionStatPlot <- function(exp.data, meta.data, stat.by, group.by = NULL, sp } if (isTRUE(same.y.lims)) { - valus <- as.matrix(dat_use[, stat.by, drop = FALSE])[is.finite(as.matrix(dat_use[, stat.by, drop = FALSE]))] + valus <- as_matrix(dat_use[, stat.by, drop = FALSE])[is.finite(as_matrix(dat_use[, stat.by, drop = FALSE]))] if (is.null(y.max)) { y.max <- max(valus, na.rm = TRUE) } else if (is.character(y.max)) { @@ -5241,7 +5243,7 @@ FeatureCorPlot <- function(srt, features, group.by = NULL, split.by = NULL, cell dat_gene <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_meta) > 0) { - dat_meta <- as.matrix(srt@meta.data[, features_meta, drop = FALSE]) + dat_meta <- as_matrix(srt@meta.data[, features_meta, drop = FALSE]) } else { dat_meta <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } @@ -5255,7 +5257,7 @@ FeatureCorPlot <- function(srt, features, group.by = NULL, split.by = NULL, cell stop("'features' must be type of numeric variable.") } if (!inherits(dat_exp, "dgCMatrix")) { - dat_exp <- as.sparse(as.matrix(dat_exp)) + dat_exp <- as.sparse(as_matrix(dat_exp)) } if (length(features) > 10 && !isTRUE(force)) { warning("More than 10 features to be paired compared which will generate more than 50 plots.", immediate. = TRUE) @@ -5664,7 +5666,7 @@ CellDensityPlot <- function(srt, features, group.by = NULL, split.by = NULL, ass dat_gene <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } if (length(features_meta) > 0) { - dat_meta <- as.matrix(srt@meta.data[, features_meta, drop = FALSE]) + dat_meta <- as_matrix(srt@meta.data[, features_meta, drop = FALSE]) } else { dat_meta <- matrix(nrow = ncol(srt@assays[[1]]), ncol = 0) } @@ -5688,10 +5690,10 @@ CellDensityPlot <- function(srt, features, group.by = NULL, split.by = NULL, ass } if (isTRUE(same.y.lims) && is.null(y.max)) { - y.max <- max(as.matrix(dat_exp[, features])[is.finite(as.matrix(dat_exp[, features]))], na.rm = TRUE) + y.max <- max(as_matrix(dat_exp[, features])[is.finite(as_matrix(dat_exp[, features]))], na.rm = TRUE) } if (isTRUE(same.y.lims) && is.null(y.min)) { - y.min <- min(as.matrix(dat_exp[, features])[is.finite(as.matrix(dat_exp[, features]))], na.rm = TRUE) + y.min <- min(as_matrix(dat_exp[, features])[is.finite(as_matrix(dat_exp[, features]))], na.rm = TRUE) } plist <- list() @@ -6119,7 +6121,7 @@ PAGAPlot <- function(srt, paga = srt@misc$paga, type = "connectivities", } out <- GraphPlot( - node = dat, edge = as.matrix(connectivities), node_coord = paste0(reduction_key, dims), + node = dat, edge = as_matrix(connectivities), node_coord = paste0(reduction_key, dims), node_group = groups, node_palette = node_palette, node_palcolor = node_palcolor, node_size = node_size, node_alpha = node_alpha, node_highlight = node_highlight, node_highlight_color = node_highlight_color, label = label, label.size = label.size, label.fg = label.fg, label.bg = label.bg, label.bg.r = label.bg.r, @@ -6391,7 +6393,7 @@ GraphPlot <- function(node, edge, transition = NULL, } if (!is.null(transition)) { - trans2 <- trans1 <- as.matrix(transition) + trans2 <- trans1 <- as_matrix(transition) trans1[lower.tri(trans1)] <- 0 trans2[upper.tri(trans2)] <- 0 trans <- t(trans1) - trans2 @@ -6910,7 +6912,7 @@ compute_velocity_on_grid <- function(X_emb, V_emb, gr <- seq(m, M, length.out = ceiling(50 * density)) grs <- c(grs, list(gr)) } - X_grid <- as.matrix(expand.grid(grs)) + X_grid <- as_matrix(expand.grid(grs)) d <- dist( x = as.sparse(X_emb), @@ -6918,8 +6920,8 @@ compute_velocity_on_grid <- function(X_emb, V_emb, method = "euclidean", use_nan = TRUE ) - neighbors <- t(as.matrix(apply(d, 2, function(x) order(x, decreasing = FALSE)[1:n_neighbors]))) - dists <- t(as.matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1:n_neighbors]]))) + neighbors <- t(as_matrix(apply(d, 2, function(x) order(x, decreasing = FALSE)[1:n_neighbors]))) + dists <- t(as_matrix(apply(d, 2, function(x) x[order(x, decreasing = FALSE)[1:n_neighbors]]))) # ggplot() + # annotate(geom = "point", x = X_grid[, 1], y = X_grid[, 2]) + @@ -7753,8 +7755,7 @@ mestimate <- function(data) { #' "Rbp4", "Pyy", # Endocrine #' "Ins1", "Gcg", "Sst", "Ghrl" # Beta, Alpha, Delta, Epsilon #' ), -#' group.by = c("CellType", "SubCellType"), -#' show_row_names = TRUE +#' group.by = c("CellType", "SubCellType") #' ) #' ht1$plot #' panel_fix(ht1$plot, height = 4, width = 6, raster = TRUE, dpi = 50) @@ -7765,14 +7766,12 @@ mestimate <- function(data) { #' ht2 <- GroupHeatmap( #' srt = pancreas_sub, features = de_filter$gene, group.by = "CellType", #' split.by = "Phase", cell_split_palette = "Dark2", -#' cluster_rows = TRUE, cluster_columns = TRUE, -#' nlabel = 10, show_row_names = FALSE +#' cluster_rows = TRUE, cluster_columns = TRUE #' ) #' ht2$plot #' #' ht3 <- GroupHeatmap( #' srt = pancreas_sub, features = de_filter$gene, feature_split = de_filter$group1, group.by = "CellType", -#' nlabel = 20, show_row_names = FALSE, #' species = "Mus_musculus", db = "GO_BP", anno_terms = TRUE, anno_keys = TRUE, anno_features = TRUE #' ) #' ht3$plot @@ -7787,10 +7786,10 @@ mestimate <- function(data) { #' features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", #' heatmap_palette = "YlOrRd", #' cell_annotation = c("Phase", "G2M_score", "Neurod2"), cell_annotation_palette = c("Dark2", "Paired", "Paired"), -#' cell_annotation_params = list(height = grid::unit(0.5, "in")), +#' cell_annotation_params = list(height = grid::unit(20, "mm")), #' feature_annotation = c("TF", "CSPA"), #' feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), -#' add_dot = TRUE, add_bg = TRUE, show_row_names = TRUE +#' add_dot = TRUE, add_bg = TRUE, nlabel = 0, show_row_names = TRUE #' ) #' ht4$plot #' @@ -7798,25 +7797,25 @@ mestimate <- function(data) { #' features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", #' heatmap_palette = "YlOrRd", #' cell_annotation = c("Phase", "G2M_score", "Neurod2"), cell_annotation_palette = c("Dark2", "Paired", "Paired"), -#' cell_annotation_params = list(width = grid::unit(0.5, "in")), +#' cell_annotation_params = list(width = grid::unit(20, "mm")), #' feature_annotation = c("TF", "CSPA"), #' feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), #' add_dot = TRUE, add_bg = TRUE, -#' flip = TRUE, column_title_rot = 45, show_column_names = TRUE +#' flip = TRUE, column_title_rot = 45, nlabel = 0, show_row_names = TRUE #' ) #' ht5$plot #' #' ht6 <- GroupHeatmap(pancreas_sub, #' features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", #' add_violin = TRUE, cluster_rows = TRUE, -#' show_row_names = TRUE +#' nlabel = 0, show_row_names = TRUE #' ) #' ht6$plot #' #' ht7 <- GroupHeatmap(pancreas_sub, #' features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", #' add_violin = TRUE, fill.by = "expression", fill_palette = "Blues", cluster_rows = TRUE, -#' show_row_names = TRUE +#' nlabel = 0, show_row_names = TRUE #' ) #' ht7$plot #' @@ -7824,7 +7823,8 @@ mestimate <- function(data) { #' features = de_top$gene, group.by = "CellType", split.by = "Phase", n_split = 4, #' cluster_rows = TRUE, cluster_columns = TRUE, cluster_row_slices = TRUE, cluster_column_slices = TRUE, #' add_dot = TRUE, add_reticle = TRUE, heatmap_palette = "viridis", -#' show_row_names = TRUE, ht_params = list(row_gap = grid::unit(0, "mm"), row_names_gp = grid::gpar(fontsize = 10)) +#' nlabel = 0, show_row_names = TRUE, +#' ht_params = list(row_gap = grid::unit(0, "mm"), row_names_gp = grid::gpar(fontsize = 10)) #' ) #' ht8$plot #' @@ -7862,8 +7862,8 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, add_violin = FALSE, fill.by = "feature", fill_palette = "Dark2", fill_palcolor = NULL, heatmap_palette = "RdBu", heatmap_palcolor = NULL, group_palette = "Paired", group_palcolor = NULL, cell_split_palette = "simspec", cell_split_palcolor = NULL, feature_split_palette = "simspec", feature_split_palcolor = NULL, - cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = if (flip) list(width = grid::unit(1, "cm")) else list(height = grid::unit(1, "cm")), - feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = list(), + cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else list(height = grid::unit(20, "mm")), + feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else list(width = grid::unit(5, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, height = NULL, width = NULL, units = "inch", seed = 11, ht_params = list()) { set.seed(seed) @@ -8083,7 +8083,7 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, gene_unique <- features_unique[features %in% rownames(srt@assays[[assay]])] meta <- features[features %in% colnames(srt@meta.data)] - mat_raw <- as.matrix(rbind(slot(srt@assays[[assay]], slot)[gene, cells, drop = FALSE], t(srt@meta.data[cells, meta, drop = FALSE])))[features, , drop = FALSE] + mat_raw <- as_matrix(rbind(slot(srt@assays[[assay]], slot)[gene, cells, drop = FALSE], t(srt@meta.data[cells, meta, drop = FALSE])))[features, , drop = FALSE] rownames(mat_raw) <- features_unique if (isTRUE(lib_normalize) && min(mat_raw, na.rm = TRUE) >= 0) { if (!is.null(libsize)) { @@ -8309,8 +8309,8 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, subplots <- CellStatPlot(srt, flip = flip, cells = names(cell_groups[[cell_group]]), plot_type = "pie", - group.by = cell_group, stat.by = cellan, split.by = split.by, - palette = palette, palcolor = palcolor, + stat.by = cellan, group.by = cell_group, split.by = split.by, + palette = palette, palcolor = palcolor, title = NULL, individual = TRUE, combine = FALSE ) subplots_list[[paste0(cellan, ":", cell_group)]] <- subplots @@ -8318,7 +8318,7 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", cell_group, "-", nm), "'; grid.draw(g) " @@ -8374,7 +8374,7 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(legend.position = 'none')); + g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", cell_group, "-", nm), "'; grid.draw(g) " @@ -8710,7 +8710,7 @@ GroupHeatmap <- function(srt, features = NULL, group.by = NULL, split.by = NULL, # legend <- get_legend(vlnplots[[1]]) # funbody <- paste0( # " - # g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + # g <- as_grob(subplots_list[['", cellan, ":", cell_group, "']]", "[['", nm, "']] + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); # grid.draw(g) # " # ) @@ -9063,8 +9063,8 @@ FeatureHeatmap <- function(srt, features = NULL, cells = NULL, group.by = NULL, nlabel = 20, features_label = NULL, label_size = 10, label_color = "black", heatmap_palette = "RdBu", heatmap_palcolor = NULL, group_palette = "Paired", group_palcolor = NULL, cell_split_palette = "simspec", cell_split_palcolor = NULL, feature_split_palette = "simspec", feature_split_palcolor = NULL, - cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = list(), - feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = list(), + cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = if (flip) list(width = grid::unit(5, "mm")) else list(height = grid::unit(5, "mm")), + feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else list(width = grid::unit(5, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, height = NULL, width = NULL, units = "inch", seed = 11, ht_params = list()) { set.seed(seed) @@ -9270,7 +9270,7 @@ FeatureHeatmap <- function(srt, features = NULL, cells = NULL, group.by = NULL, gene_unique <- features_unique[features %in% rownames(srt@assays[[assay]])] meta <- features[features %in% colnames(srt@meta.data)] all_cells <- unique(unlist(lapply(cell_groups, names))) - mat_raw <- as.matrix(rbind(slot(srt@assays[[assay]], slot)[gene, all_cells, drop = FALSE], t(srt@meta.data[all_cells, meta, drop = FALSE])))[features, , drop = FALSE] + mat_raw <- as_matrix(rbind(slot(srt@assays[[assay]], slot)[gene, all_cells, drop = FALSE], t(srt@meta.data[all_cells, meta, drop = FALSE])))[features, , drop = FALSE] rownames(mat_raw) <- features_unique if (isTRUE(lib_normalize) && min(mat_raw, na.rm = TRUE) >= 0) { if (!is.null(libsize)) { @@ -10016,7 +10016,7 @@ FeatureCorHeatmap <- function(srt, features, cells) { #' query_group = "SubCellType", ref_group = "celltype", #' query_cell_annotation = "Phase", query_cell_annotation_palette = "Set2", #' ref_cell_annotation = "tech", ref_cell_annotation_palette = "Set3", -#' width = 3, height = 2 +#' width = 4, height = 4 #' ) #' ht2$plot #' @@ -10064,8 +10064,8 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, heatmap_palette = "RdBu", heatmap_palcolor = NULL, query_group_palette = "Paired", query_group_palcolor = NULL, ref_group_palette = "simspec", ref_group_palcolor = NULL, - query_cell_annotation = NULL, query_cell_annotation_palette = "Paired", query_cell_annotation_palcolor = NULL, query_cell_annotation_params = if (flip) list(height = grid::unit(1, "cm")) else list(width = grid::unit(1, "cm")), - ref_cell_annotation = NULL, ref_cell_annotation_palette = "Paired", ref_cell_annotation_palcolor = NULL, ref_cell_annotation_params = if (flip) list(width = grid::unit(1, "cm")) else list(height = grid::unit(1, "cm")), + query_cell_annotation = NULL, query_cell_annotation_palette = "Paired", query_cell_annotation_palcolor = NULL, query_cell_annotation_params = if (flip) list(height = grid::unit(20, "mm")) else list(width = grid::unit(20, "mm")), + ref_cell_annotation = NULL, ref_cell_annotation_palette = "Paired", ref_cell_annotation_palcolor = NULL, ref_cell_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else list(height = grid::unit(20, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, height = NULL, width = NULL, units = "inch", seed = 11, ht_params = list()) { set.seed(seed) @@ -10127,10 +10127,10 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, distance_matrix <- srt_query@tools[[paste0(prefix, "_classification")]][["distance_matrix"]] distance_metric <- srt_query@tools[[paste0(prefix, "_classification")]][["distance_metric"]] if (distance_metric %in% simil_method) { - simil_matrix <- t(as.matrix(1 - distance_matrix)) + simil_matrix <- t(as_matrix(1 - distance_matrix)) simil_name <- paste0(capitalize(distance_metric), " similarity") } else if (distance_metric %in% dist_method) { - simil_matrix <- t(as.matrix(1 - distance_matrix / max(distance_matrix, na.rm = TRUE))) + simil_matrix <- t(as_matrix(1 - distance_matrix / max(distance_matrix, na.rm = TRUE))) simil_name <- paste0("1-dist[", distance_metric, "]/max(dist[", distance_metric, "])") } simil_matrix[is.infinite(simil_matrix)] <- max(abs(simil_matrix[!is.infinite(simil_matrix)]), na.rm = TRUE) * ifelse(simil_matrix[is.infinite(simil_matrix)] > 0, 1, -1) @@ -10345,7 +10345,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(query_subplots_list[['", cellan, ":", query_group, "']]", "[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + g <- as_grob(query_subplots_list[['", cellan, ":", query_group, "']]", "[['", nm, "']] + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", query_group, "-", nm), "'; grid.draw(g) " @@ -10423,7 +10423,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(query_subplots_list[['", cellan, ":", query_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(legend.position = 'none')); + g <- as_grob(query_subplots_list[['", cellan, ":", query_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", query_group, "-", nm), "'; grid.draw(g) " @@ -10518,7 +10518,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(ref_subplots_list[['", cellan, ":", ref_group, "']]", "[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + g <- as_grob(ref_subplots_list[['", cellan, ":", ref_group, "']]", "[['", nm, "']] + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", ref_group, "-", nm), "'; grid.draw(g) " @@ -10596,7 +10596,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, for (nm in names(subplots)) { funbody <- paste0( " - g <- as_grob(ref_subplots_list[['", cellan, ":", ref_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(legend.position = 'none')); + g <- as_grob(ref_subplots_list[['", cellan, ":", ref_group, "']]", "[['", nm, "']] + facet_null() + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", paste0(cellan, ":", ref_group, "-", nm), "'; grid.draw(g) " @@ -10892,7 +10892,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, #' cell_annotation_palette = c("Paired", "simspec", "Purples"), #' separate_annotation = list("SubCellType", c("Nnat", "Irx1")), #' separate_annotation_palette = c("Paired", "Set1"), -#' separate_annotation_params = list(height = grid::unit(2, "cm")), +#' separate_annotation_params = list(height = grid::unit(20, "mm")), #' feature_annotation = c("TF", "CSPA"), #' feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), #' pseudotime_label = 25, @@ -10912,7 +10912,7 @@ CellCorHeatmap <- function(srt_query, srt_ref = NULL, bulk_ref = NULL, #' cell_annotation_palette = c("Paired", "simspec", "Purples"), #' separate_annotation = list("SubCellType", c("Nnat", "Irx1")), #' separate_annotation_palette = c("Paired", "Set1"), -#' separate_annotation_params = list(width = grid::unit(2, "cm")), +#' separate_annotation_params = list(width = grid::unit(20, "mm")), #' feature_annotation = c("TF", "CSPA"), #' feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), #' pseudotime_label = 25, @@ -10967,9 +10967,9 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b heatmap_palette = "viridis", heatmap_palcolor = NULL, pseudotime_palette = "cividis", pseudotime_palcolor = NULL, feature_split_palette = "simspec", feature_split_palcolor = NULL, - cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = list(), - feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = list(), - separate_annotation = NULL, separate_annotation_palette = "Paired", separate_annotation_palcolor = NULL, separate_annotation_params = if (flip) list(width = grid::unit(2, "cm")) else list(height = grid::unit(2, "cm")), + cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, cell_annotation_params = if (flip) list(width = grid::unit(5, "mm")) else list(height = grid::unit(5, "mm")), + feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else list(width = grid::unit(5, "mm")), + separate_annotation = NULL, separate_annotation_palette = "Paired", separate_annotation_palcolor = NULL, separate_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else list(height = grid::unit(20, "mm")), reverse_ht = NULL, use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, height = NULL, width = NULL, units = "inch", seed = 11, ht_params = list()) { set.seed(seed) @@ -11208,7 +11208,7 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b Y_libsize <- colSums(slot(srt@assays[[assay]], "counts")) for (l in lineages) { cells <- gsub(pattern = l, replacement = "", x = cell_order_list[[l]]) - mat_tmp <- as.matrix(rbind(slot(srt@assays[[assay]], slot)[gene, cells, drop = FALSE], t(srt@meta.data[cells, meta, drop = FALSE])))[features, , drop = FALSE] + mat_tmp <- as_matrix(rbind(slot(srt@assays[[assay]], slot)[gene, cells, drop = FALSE], t(srt@meta.data[cells, meta, drop = FALSE])))[features, , drop = FALSE] if (isTRUE(lib_normalize) && min(mat_tmp, na.rm = TRUE) >= 0) { if (!is.null(libsize)) { libsize_use <- libsize @@ -11375,7 +11375,7 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b nm <- paste0(cellan, ":", l) funbody <- paste0( " - g <- as_grob(subplots_list[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + g <- as_grob(subplots_list[['", nm, "']] + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", nm, "'; grid.draw(g); grid.rect(gp = gpar(fill = 'transparent', col = 'black')); @@ -11421,7 +11421,7 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b nm <- paste0(paste0(cellan, collapse = ","), ":", l) funbody <- paste0( " - g <- as_grob(subplots_list[['", nm, "']] + theme_void() + theme(legend.position = 'none')); + g <- as_grob(subplots_list[['", nm, "']] + theme_void() + theme(plot.title = element_blank(), plot.subtitle = element_blank(), legend.position = 'none')); g$name <- '", nm, "'; grid.draw(g); grid.rect(gp = gpar(fill = 'transparent', col = 'black')); @@ -11650,8 +11650,8 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b ha_list[[l]] <- anno_simple( x = is.na(feature_metadata[, paste0(l, order_by)]) + 0, col = c("0" = "#181830", "1" = "transparent"), - width = unit(0.5, "cm"), - height = unit(0.5, "cm"), + width = unit(5, "mm"), + height = unit(5, "mm"), which = ifelse(flip, "column", "row") ) } @@ -11802,9 +11802,7 @@ DynamicHeatmap <- function(srt, lineages, features = NULL, use_fitted = FALSE, b if ((!is.null(row_split) && length(index) > 0) || any(c(anno_terms, anno_keys, anno_features)) || !is.null(width) || !is.null(height)) { fix <- TRUE if (is.null(width) || is.null(height)) { - message("The size of the heatmap is fixed because certain elements are not scalable.\n - The width and height of the heatmap are determined by the size of the current viewport.\n - If you want to have more control over the size, you can manually set the parameters 'width' and 'height'.") + message("\nThe size of the heatmap is fixed because certain elements are not scalable.\nThe width and height of the heatmap are determined by the size of the current viewport.\nIf you want to have more control over the size, you can manually set the parameters 'width' and 'height'.\n") } } else { fix <- FALSE @@ -12088,10 +12086,10 @@ DynamicPlot <- function(srt, lineages, features, group.by = NULL, cells = NULL, upr_matrix <- cbind(upr_matrix, srt_tmp@tools[[paste0("DynamicFeatures_", l)]][["upr_matrix"]][, feature_calcu, drop = FALSE]) lwr_matrix <- cbind(lwr_matrix, srt_tmp@tools[[paste0("DynamicFeatures_", l)]][["lwr_matrix"]][, feature_calcu, drop = FALSE]) } - raw_matrix_list[[l]] <- as.matrix(raw_matrix[, features, drop = FALSE]) - fitted_matrix_list[[l]] <- as.matrix(fitted_matrix[, features, drop = FALSE]) - upr_matrix_list[[l]] <- as.matrix(upr_matrix[, features, drop = FALSE]) - lwr_matrix_list[[l]] <- as.matrix(lwr_matrix[, features, drop = FALSE]) + raw_matrix_list[[l]] <- as_matrix(raw_matrix[, features, drop = FALSE]) + fitted_matrix_list[[l]] <- as_matrix(fitted_matrix[, features, drop = FALSE]) + upr_matrix_list[[l]] <- as_matrix(upr_matrix[, features, drop = FALSE]) + lwr_matrix_list[[l]] <- as_matrix(lwr_matrix[, features, drop = FALSE]) cell_union <- unique(c(cell_union, rownames(raw_matrix))) } @@ -13394,7 +13392,7 @@ adjustlayout <- function(graph, layout, width, height = 2, scale = 100, iter = 1 # } for (i in seq_len(iter)) { - dist_matrix <- as.matrix(dist(layout)) + dist_matrix <- as_matrix(dist(layout)) nearest_neighbors <- apply(dist_matrix, 2, function(x) which(x == min(x[x > 0])), simplify = FALSE) # nearest_neighbors <- apply(dist_matrix, 2, function(x) { # head(order(x), 3)[-1] diff --git a/R/SCP-workflow.R b/R/SCP-workflow.R index 66165414..f6897fa2 100644 --- a/R/SCP-workflow.R +++ b/R/SCP-workflow.R @@ -1600,7 +1600,7 @@ scVI_integrate <- function(srtMerge = NULL, batch = NULL, append = TRUE, srtList model$train() srtIntegrated <- srtMerge srtMerge <- NULL - corrected <- t(as.matrix(model$get_normalized_expression())) + corrected <- t(as_matrix(model$get_normalized_expression())) srtIntegrated[["scVIcorrected"]] <- CreateAssayObject(counts = corrected) DefaultAssay(srtIntegrated) <- "scVIcorrected" VariableFeatures(srtIntegrated[["scVIcorrected"]]) <- HVF @@ -1619,7 +1619,7 @@ scVI_integrate <- function(srtMerge = NULL, batch = NULL, append = TRUE, srtList srtMerge <- NULL } - latent <- as.matrix(model$get_latent_representation()) + latent <- as_matrix(model$get_latent_representation()) rownames(latent) <- colnames(srtIntegrated) colnames(latent) <- paste0("scVI_", seq_len(ncol(latent))) srtIntegrated[["scVI"]] <- CreateDimReducObject(embeddings = latent, key = "scVI_", assay = DefaultAssay(srtIntegrated)) @@ -1990,7 +1990,7 @@ fastMNN_integrate <- function(srtMerge = NULL, batch = NULL, append = TRUE, srtL srtIntegrated <- srtMerge srtMerge <- NULL - srtIntegrated[["fastMNNcorrected"]] <- CreateAssayObject(counts = as.matrix(out@assays@data$reconstructed)) + srtIntegrated[["fastMNNcorrected"]] <- CreateAssayObject(counts = as_matrix(out@assays@data$reconstructed)) DefaultAssay(srtIntegrated) <- "fastMNNcorrected" VariableFeatures(srtIntegrated[["fastMNNcorrected"]]) <- HVF reduction <- out@int_colData$reducedDims$corrected @@ -2352,7 +2352,7 @@ Scanorama_integrate <- function(srtMerge = NULL, batch = NULL, append = TRUE, sr assaylist <- list() genelist <- list() for (i in seq_along(srtList)) { - assaylist[[i]] <- t(as.matrix(GetAssayData(object = srtList[[i]], slot = "data", assay = DefaultAssay(srtList[[i]]))[HVF, , drop = FALSE])) + assaylist[[i]] <- t(as_matrix(GetAssayData(object = srtList[[i]], slot = "data", assay = DefaultAssay(srtList[[i]]))[HVF, , drop = FALSE])) genelist[[i]] <- HVF } if (isTRUE(return_corrected)) { diff --git a/R/Seurat-function.R b/R/Seurat-function.R index a7b09cf7..460b7788 100644 --- a/R/Seurat-function.R +++ b/R/Seurat-function.R @@ -134,7 +134,7 @@ RunNMF.default <- function(object, assay = NULL, slot = "data", nbes = 50, if (nmf.method == "NMF") { check_R("NMF") seed <- NMF::seed - nmf.results <- NMF::nmf(x = as.matrix(t(object)), rank = nbes) + nmf.results <- NMF::nmf(x = as_matrix(t(object)), rank = nbes) cell.embeddings <- nmf.results@fit@W feature.loadings <- t(nmf.results@fit@H) } @@ -271,7 +271,7 @@ RunMDS.default <- function(object, assay = NULL, slot = "data", object <- t(object) } nmds <- min(nmds, nrow(x = object) - 1) - x <- t(as.matrix(object)) + x <- t(as_matrix(object)) cell.dist <- as.dist(dist(x = x, method = dist.method)) if (mds.method == "cmdscale") { mds.results <- cmdscale(cell.dist, k = nmds, eig = TRUE, ...) @@ -402,8 +402,8 @@ RunGLMPCA.default <- function(object, assay = NULL, slot = "counts", fam <- match.arg(fam) glmpca_results <- glmpca::glmpca(Y = object, L = L, fam = fam, ...) glmpca_dimnames <- paste0(reduction.key, seq_len(L)) - factors <- as.matrix(glmpca_results$factors) - loadings <- as.matrix(glmpca_results$loadings) + factors <- as_matrix(glmpca_results$factors) + loadings <- as_matrix(glmpca_results$loadings) colnames(x = factors) <- glmpca_dimnames colnames(x = loadings) <- glmpca_dimnames factors_l2_norm <- sqrt(colSums(factors^2)) @@ -471,7 +471,7 @@ RunDM.Seurat <- function(object, verbose = TRUE, seed.use = 11, ...) { if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < ndcs) { stop("Please provide as many or more features than ndcs: ", length(x = features), " features provided, ", @@ -523,7 +523,7 @@ RunDM.default <- function(object, assay = NULL, slot = "data", set.seed(seed = seed.use) } - dm.results <- destiny::DiffusionMap(data = as.matrix(object), n_eigs = ndcs, sigma = sigma, k = k, distance = dist.method, verbose = verbose, ...) + dm.results <- destiny::DiffusionMap(data = as_matrix(object), n_eigs = ndcs, sigma = sigma, k = k, distance = dist.method, verbose = verbose, ...) cell.embeddings <- dm.results@eigenvectors rownames(x = cell.embeddings) <- rownames(object) @@ -599,7 +599,7 @@ RunUMAP2.Seurat <- function(object, } if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < n.components) { stop( "Please provide as many or more features than n.components: ", @@ -611,7 +611,7 @@ RunUMAP2.Seurat <- function(object, ) } } else if (!is.null(x = dims)) { - data.use <- as.matrix(Embeddings(object[[reduction]])[, dims]) + data.use <- as_matrix(Embeddings(object[[reduction]])[, dims]) assay <- DefaultAssay(object = object[[reduction]]) if (length(x = dims) < n.components) { stop( @@ -793,7 +793,7 @@ RunUMAP2.default <- function(object, assay = NULL, } else { obs_sample <- 1:ncol(object) } - if (!isSymmetric(as.matrix(object[obs_sample, obs_sample]))) { + if (!isSymmetric(as_matrix(object[obs_sample, obs_sample]))) { stop("Graph must be a symmetric matrix.") } @@ -810,7 +810,7 @@ RunUMAP2.default <- function(object, assay = NULL, # if (inherits(object, what = "dgCMatrix")) { # object <- as_matrix(object) # } else { - # object <- as.matrix(object) + # object <- as_matrix(object) # } # if (!isSymmetric(object)) { # stop("Graph must be a symmetric matrix.") @@ -941,7 +941,7 @@ RunUMAP2.default <- function(object, assay = NULL, } else { obs_sample <- 1:ncol(object) } - if (!isSymmetric(as.matrix(object[obs_sample, obs_sample]))) { + if (!isSymmetric(as_matrix(object[obs_sample, obs_sample]))) { stop("Graph must be a symmetric matrix.") } val <- split(object@x, rep(1:ncol(object), diff(object@p))) @@ -959,8 +959,8 @@ RunUMAP2.default <- function(object, assay = NULL, }, x = val, y = val)) idx[is.na(idx)] <- sample(1:nrow(object), size = sum(is.na(idx)), replace = TRUE) nn <- list(idx = idx, dist = max(connectivity) - connectivity + min(diff(range(connectivity)), 1) / 1e50) - # idx <- t(as.matrix(apply(object, 2, function(x) order(x, decreasing = TRUE)[1:n.neighbors]))) - # connectivity <- t(as.matrix(apply(object, 2, function(x) x[order(x, decreasing = TRUE)[1:n.neighbors]]))) + # idx <- t(as_matrix(apply(object, 2, function(x) order(x, decreasing = TRUE)[1:n.neighbors]))) + # connectivity <- t(as_matrix(apply(object, 2, function(x) x[order(x, decreasing = TRUE)[1:n.neighbors]]))) out <- uwot::umap( X = NULL, nn_method = nn, n_threads = 1, n_components = n.components, metric = metric, n_epochs = n.epochs, learning_rate = learning.rate, @@ -1080,8 +1080,8 @@ RunUMAP2.default <- function(object, assay = NULL, return(reduction) } if (inherits(x = object, what = "Graph")) { - match_k <- t(as.matrix(apply(object, 2, function(x) order(x, decreasing = TRUE)[1:n.neighbors]))) - match_k_connectivity <- t(as.matrix(apply(object, 2, function(x) x[order(x, decreasing = TRUE)[1:n.neighbors]]))) + match_k <- t(as_matrix(apply(object, 2, function(x) order(x, decreasing = TRUE)[1:n.neighbors]))) + match_k_connectivity <- t(as_matrix(apply(object, 2, function(x) x[order(x, decreasing = TRUE)[1:n.neighbors]]))) object <- list(idx = match_k, dist = max(match_k_connectivity) - match_k_connectivity) if (is.null(model$num_precomputed_nns) || model$num_precomputed_nns == 0) { model$num_precomputed_nns <- 1 @@ -1168,7 +1168,7 @@ RunPaCMAP.Seurat <- function(object, reduction = "pca", dims = NULL, features = } if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < n_components) { stop( "Please provide as many or more features than n_components: ", @@ -1297,7 +1297,7 @@ RunPHATE.Seurat <- function(object, reduction = "pca", dims = NULL, features = N } if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < n_components) { stop( "Please provide as many or more features than n_components: ", @@ -1448,7 +1448,7 @@ RunTriMap.Seurat <- function(object, reduction = "pca", dims = NULL, features = } if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < n_components) { stop( "Please provide as many or more features than n_components: ", @@ -1570,7 +1570,7 @@ RunLargeVis.Seurat <- function(object, reduction = "pca", dims = NULL, features } if (!is.null(x = features)) { assay <- assay %||% DefaultAssay(object = object) - data.use <- as.matrix(x = t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) + data.use <- as_matrix(t(x = GetAssayData(object = object, slot = slot, assay = assay)[features, , drop = FALSE])) if (ncol(x = data.use) < n_components) { stop( "Please provide as many or more features than n_components: ", @@ -1824,7 +1824,7 @@ RunHarmony2.Seurat <- function(object, group.by.vars, ... ) - harmonyEmbed <- t(as.matrix(harmonyObject$Z_corr)) + harmonyEmbed <- t(as_matrix(harmonyObject$Z_corr)) rownames(harmonyEmbed) <- row.names(data.use) colnames(harmonyEmbed) <- paste0(reduction.name, "_", seq_len(ncol(harmonyEmbed))) diff --git a/R/utils.R b/R/utils.R index 567b01f0..f0c043f2 100644 --- a/R/utils.R +++ b/R/utils.R @@ -816,24 +816,27 @@ unnest <- function(data, cols, keep_empty = FALSE) { } #' Attempts to turn a dgCMatrix into a dense matrix -#' @param matrix A dgCMatrix +#' +#' @examples +#' data("pancreas_sub") +#' system.time(mat1 <- as.matrix(slot(pancreas_sub[["RNA"]], "counts"))) +#' system.time(mat2 <- as_matrix(slot(pancreas_sub[["RNA"]], "counts"))) +#' identical(mat1, mat2) +#' +#' @param x A matrix. #' @useDynLib SCP +#' @importFrom Matrix as.matrix #' @export -as_matrix <- function(matrix) { +as_matrix <- function(x) { if (!inherits(matrix, "dgCMatrix")) { - stop("matrix is not a dgCMatrix.") + return(as.matrix(x)) + } else { + row_pos <- x@i + col_pos <- findInterval(seq_along(x@x) - 1, x@p[-1]) + out <- asMatrix(rp = row_pos, cp = col_pos, z = x@x, nrows = x@Dim[1], ncols = x@Dim[2]) + attr(out, "dimnames") <- list(x@Dimnames[[1]], x@Dimnames[[2]]) + return(out) } - row_pos <- matrix@i - col_pos <- findInterval(seq(matrix@x) - 1, matrix@p[-1]) - - out <- asMatrix( - rp = row_pos, cp = col_pos, z = matrix@x, - nrows = matrix@Dim[1], ncols = matrix@Dim[2] - ) - - row.names(out) <- matrix@Dimnames[[1]] - colnames(out) <- matrix@Dimnames[[2]] - return(out) } #' Capitalizes the characters diff --git a/man/CellCorHeatmap.Rd b/man/CellCorHeatmap.Rd index b07c141c..9c12ba94 100644 --- a/man/CellCorHeatmap.Rd +++ b/man/CellCorHeatmap.Rd @@ -59,13 +59,13 @@ CellCorHeatmap( query_cell_annotation = NULL, query_cell_annotation_palette = "Paired", query_cell_annotation_palcolor = NULL, - query_cell_annotation_params = if (flip) list(height = grid::unit(1, "cm")) else - list(width = grid::unit(1, "cm")), + query_cell_annotation_params = if (flip) list(height = grid::unit(20, "mm")) else + list(width = grid::unit(20, "mm")), ref_cell_annotation = NULL, ref_cell_annotation_palette = "Paired", ref_cell_annotation_palcolor = NULL, - ref_cell_annotation_params = if (flip) list(width = grid::unit(1, "cm")) else - list(height = grid::unit(1, "cm")), + ref_cell_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else + list(height = grid::unit(20, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, @@ -241,7 +241,7 @@ ht2 <- CellCorHeatmap( query_group = "SubCellType", ref_group = "celltype", query_cell_annotation = "Phase", query_cell_annotation_palette = "Set2", ref_cell_annotation = "tech", ref_cell_annotation_palette = "Set3", - width = 3, height = 2 + width = 4, height = 4 ) ht2$plot diff --git a/man/CellScoring.Rd b/man/CellScoring.Rd index 5d5154de..bc53574b 100644 --- a/man/CellScoring.Rd +++ b/man/CellScoring.Rd @@ -86,7 +86,7 @@ data("pancreas_sub") ccgenes <- CC_GenePrefetch("Mus_musculus") pancreas_sub <- CellScoring( srt = pancreas_sub, - features = list(S = ccgenes$cc_S_genes, G2M = ccgenes$cc_G2M_genes), + features = list(S = ccgenes$S, G2M = ccgenes$G2M), method = "Seurat", name = "CC" ) CellDimPlot(pancreas_sub, "CC_classification") diff --git a/man/DynamicHeatmap.Rd b/man/DynamicHeatmap.Rd index b929852b..c3a248c9 100644 --- a/man/DynamicHeatmap.Rd +++ b/man/DynamicHeatmap.Rd @@ -100,16 +100,18 @@ DynamicHeatmap( cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, - cell_annotation_params = list(), + cell_annotation_params = if (flip) list(width = grid::unit(5, "mm")) else list(height = + grid::unit(5, "mm")), feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, - feature_annotation_params = list(), + feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else + list(width = grid::unit(5, "mm")), separate_annotation = NULL, separate_annotation_palette = "Paired", separate_annotation_palcolor = NULL, - separate_annotation_params = if (flip) list(width = grid::unit(2, "cm")) else - list(height = grid::unit(2, "cm")), + separate_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else + list(height = grid::unit(20, "mm")), reverse_ht = NULL, use_raster = NULL, raster_device = "png", @@ -397,7 +399,7 @@ ht4 <- DynamicHeatmap( cell_annotation_palette = c("Paired", "simspec", "Purples"), separate_annotation = list("SubCellType", c("Nnat", "Irx1")), separate_annotation_palette = c("Paired", "Set1"), - separate_annotation_params = list(height = grid::unit(2, "cm")), + separate_annotation_params = list(height = grid::unit(20, "mm")), feature_annotation = c("TF", "CSPA"), feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), pseudotime_label = 25, @@ -417,7 +419,7 @@ ht5 <- DynamicHeatmap( cell_annotation_palette = c("Paired", "simspec", "Purples"), separate_annotation = list("SubCellType", c("Nnat", "Irx1")), separate_annotation_palette = c("Paired", "Set1"), - separate_annotation_params = list(width = grid::unit(2, "cm")), + separate_annotation_params = list(width = grid::unit(20, "mm")), feature_annotation = c("TF", "CSPA"), feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), pseudotime_label = 25, diff --git a/man/FeatureHeatmap.Rd b/man/FeatureHeatmap.Rd index fa5d7063..5f03cab0 100644 --- a/man/FeatureHeatmap.Rd +++ b/man/FeatureHeatmap.Rd @@ -93,11 +93,13 @@ FeatureHeatmap( cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, - cell_annotation_params = list(), + cell_annotation_params = if (flip) list(width = grid::unit(5, "mm")) else list(height = + grid::unit(5, "mm")), feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, - feature_annotation_params = list(), + feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else + list(width = grid::unit(5, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, diff --git a/man/FeatureStatPlot.Rd b/man/FeatureStatPlot.Rd index fbb0a04d..7697c3f3 100644 --- a/man/FeatureStatPlot.Rd +++ b/man/FeatureStatPlot.Rd @@ -32,7 +32,7 @@ FeatureStatPlot( pt.color = "grey30", pt.size = NULL, pt.alpha = 1, - jitter.width = 0.5, + jitter.width = 0.4, jitter.height = 0.1, add_trend = FALSE, trend_color = "black", @@ -266,6 +266,7 @@ FeatureStatPlot(pancreas_sub, group.by = "SubCellType", bg.by = "CellType", stack = TRUE, flip = TRUE ) \%>\% panel_fix_overall(width = 8, height = 5) # As the plot is created by combining, we can adjust the overall height and width directly. +FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "group") FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "feature") FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", plot.by = "feature", @@ -286,7 +287,7 @@ FeatureStatPlot(pancreas_sub, stat.by = c( library(Matrix) data <- pancreas_sub@assays$RNA@data -pancreas_sub@assays$RNA@scale.data <- as.matrix(data / rowMeans(data)) +pancreas_sub@assays$RNA@scale.data <- as_matrix(data / rowMeans(data)) FeatureStatPlot(pancreas_sub, stat.by = c("Neurog3", "Rbp4", "Ins1"), group.by = "CellType", slot = "scale.data", ylab = "FoldChange", same.y.lims = TRUE, y.max = 4 diff --git a/man/FetchH5.Rd b/man/FetchH5.Rd index 439970dc..9ee51669 100644 --- a/man/FetchH5.Rd +++ b/man/FetchH5.Rd @@ -43,7 +43,13 @@ This function fetches data from an hdf5 file. It can fetch gene expression data, data("pancreas_sub") pancreas_sub <- Standard_SCP(pancreas_sub) PrepareSCExplorer(pancreas_sub, base_dir = "./SCExplorer") -srt <- FetchH5(DataFile = "./SCExplorer/Data.hdf5", MetaFile = "./SCExplorer/Meta.hdf5", features = c("Ins1", "Ghrl"), metanames = c("SubCellType", "Phase"), reduction = "UMAP") +srt <- FetchH5( + DataFile = "./SCExplorer/Data.hdf5", + MetaFile = "./SCExplorer/Meta.hdf5", + features = c("Ins1", "Ghrl"), + metanames = c("SubCellType", "Phase"), + reduction = "UMAP" +) CellDimPlot(srt, group.by = c("SubCellType", "Phase"), reduction = "UMAP") FeatureDimPlot(srt, features = c("Ins1", "Ghrl"), reduction = "UMAP") } diff --git a/man/GeneConvert.Rd b/man/GeneConvert.Rd index b9a4fec3..b18cb020 100644 --- a/man/GeneConvert.Rd +++ b/man/GeneConvert.Rd @@ -87,7 +87,7 @@ homologs_counts <- aggregate( by = list(res$geneID_expand[, "symbol"]), FUN = sum ) rownames(homologs_counts) <- homologs_counts[, 1] -homologs_counts <- as(as.matrix(homologs_counts[, -1]), "dgCMatrix") +homologs_counts <- as(as_matrix(homologs_counts[, -1]), "dgCMatrix") homologs_counts } diff --git a/man/GroupHeatmap.Rd b/man/GroupHeatmap.Rd index d454313f..25e57f65 100644 --- a/man/GroupHeatmap.Rd +++ b/man/GroupHeatmap.Rd @@ -105,12 +105,13 @@ GroupHeatmap( cell_annotation = NULL, cell_annotation_palette = "Paired", cell_annotation_palcolor = NULL, - cell_annotation_params = if (flip) list(width = grid::unit(1, "cm")) else list(height = - grid::unit(1, "cm")), + cell_annotation_params = if (flip) list(width = grid::unit(20, "mm")) else list(height + = grid::unit(20, "mm")), feature_annotation = NULL, feature_annotation_palette = "Dark2", feature_annotation_palcolor = NULL, - feature_annotation_params = list(), + feature_annotation_params = if (flip) list(height = grid::unit(5, "mm")) else + list(width = grid::unit(5, "mm")), use_raster = NULL, raster_device = "png", raster_by_magick = FALSE, @@ -373,8 +374,7 @@ ht1 <- GroupHeatmap(pancreas_sub, "Rbp4", "Pyy", # Endocrine "Ins1", "Gcg", "Sst", "Ghrl" # Beta, Alpha, Delta, Epsilon ), - group.by = c("CellType", "SubCellType"), - show_row_names = TRUE + group.by = c("CellType", "SubCellType") ) ht1$plot panel_fix(ht1$plot, height = 4, width = 6, raster = TRUE, dpi = 50) @@ -385,14 +385,12 @@ de_filter <- filter(pancreas_sub@tools$DEtest_CellType$AllMarkers_wilcox, p_val_ ht2 <- GroupHeatmap( srt = pancreas_sub, features = de_filter$gene, group.by = "CellType", split.by = "Phase", cell_split_palette = "Dark2", - cluster_rows = TRUE, cluster_columns = TRUE, - nlabel = 10, show_row_names = FALSE + cluster_rows = TRUE, cluster_columns = TRUE ) ht2$plot ht3 <- GroupHeatmap( srt = pancreas_sub, features = de_filter$gene, feature_split = de_filter$group1, group.by = "CellType", - nlabel = 20, show_row_names = FALSE, species = "Mus_musculus", db = "GO_BP", anno_terms = TRUE, anno_keys = TRUE, anno_features = TRUE ) ht3$plot @@ -407,10 +405,10 @@ ht4 <- GroupHeatmap(pancreas_sub, features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", heatmap_palette = "YlOrRd", cell_annotation = c("Phase", "G2M_score", "Neurod2"), cell_annotation_palette = c("Dark2", "Paired", "Paired"), - cell_annotation_params = list(height = grid::unit(0.5, "in")), + cell_annotation_params = list(height = grid::unit(20, "mm")), feature_annotation = c("TF", "CSPA"), feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), - add_dot = TRUE, add_bg = TRUE, show_row_names = TRUE + add_dot = TRUE, add_bg = TRUE, nlabel = 0, show_row_names = TRUE ) ht4$plot @@ -418,25 +416,25 @@ ht5 <- GroupHeatmap(pancreas_sub, features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", heatmap_palette = "YlOrRd", cell_annotation = c("Phase", "G2M_score", "Neurod2"), cell_annotation_palette = c("Dark2", "Paired", "Paired"), - cell_annotation_params = list(width = grid::unit(0.5, "in")), + cell_annotation_params = list(width = grid::unit(20, "mm")), feature_annotation = c("TF", "CSPA"), feature_annotation_palcolor = list(c("gold", "steelblue"), c("forestgreen")), add_dot = TRUE, add_bg = TRUE, - flip = TRUE, column_title_rot = 45, show_column_names = TRUE + flip = TRUE, column_title_rot = 45, nlabel = 0, show_row_names = TRUE ) ht5$plot ht6 <- GroupHeatmap(pancreas_sub, features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", add_violin = TRUE, cluster_rows = TRUE, - show_row_names = TRUE + nlabel = 0, show_row_names = TRUE ) ht6$plot ht7 <- GroupHeatmap(pancreas_sub, features = de_top$gene, feature_split = de_top$group1, group.by = "CellType", add_violin = TRUE, fill.by = "expression", fill_palette = "Blues", cluster_rows = TRUE, - show_row_names = TRUE + nlabel = 0, show_row_names = TRUE ) ht7$plot @@ -444,7 +442,8 @@ ht8 <- GroupHeatmap(pancreas_sub, features = de_top$gene, group.by = "CellType", split.by = "Phase", n_split = 4, cluster_rows = TRUE, cluster_columns = TRUE, cluster_row_slices = TRUE, cluster_column_slices = TRUE, add_dot = TRUE, add_reticle = TRUE, heatmap_palette = "viridis", - show_row_names = TRUE, ht_params = list(row_gap = grid::unit(0, "mm"), row_names_gp = grid::gpar(fontsize = 10)) + nlabel = 0, show_row_names = TRUE, + ht_params = list(row_gap = grid::unit(0, "mm"), row_names_gp = grid::gpar(fontsize = 10)) ) ht8$plot diff --git a/man/RunSCExplorer.Rd b/man/RunSCExplorer.Rd index 03ca1990..3af7b87a 100644 --- a/man/RunSCExplorer.Rd +++ b/man/RunSCExplorer.Rd @@ -15,14 +15,14 @@ RunSCExplorer( initial_feature = NULL, initial_assay = NULL, initial_slot = NULL, - initial_label = "No", + initial_label = FALSE, initial_cell_palette = "Paired", initial_feature_palette = "Spectral", initial_theme = "theme_scp", initial_size = 4, initial_ncol = 3, - initial_arrange = "Row", - initial_raster = "No", + initial_arrange = NULL, + initial_raster = NULL, session_workers = 2, plotting_workers = 8, create_script = TRUE, @@ -52,7 +52,7 @@ RunSCExplorer( \item{initial_slot}{A string. The initial slot to be loaded into the app. Default is NULL.} -\item{initial_label}{A string. Whether to add labels in the initial plot. Default is "No".} +\item{initial_label}{A string. Whether to add labels in the initial plot. Default is FALSE.} \item{initial_cell_palette}{A string. The initial color palette for cells. Default is "Paired".} @@ -64,7 +64,9 @@ RunSCExplorer( \item{initial_ncol}{A numeric. The initial number of columns for arranging plots. Default is 3.} -\item{initial_arrange}{A string. The initial arrangement of plots. Default is "Row".} +\item{initial_arrange}{A logical. Whether to use "Row" as the initial arrangement. Default is TRUE.} + +\item{initial_raster}{A logical. Whether to perform rasterization in the initial plot. By default, it is set to automatic, meaning it will be TRUE if the number of cells in the initial datasets exceeds 100,000.} \item{session_workers}{A numeric. The number of workers for concurrent execution in an asynchronous programming session. Default is 2.} diff --git a/man/as_matrix.Rd b/man/as_matrix.Rd index 1687aedd..c284c0a9 100644 --- a/man/as_matrix.Rd +++ b/man/as_matrix.Rd @@ -4,11 +4,18 @@ \alias{as_matrix} \title{Attempts to turn a dgCMatrix into a dense matrix} \usage{ -as_matrix(matrix) +as_matrix(x) } \arguments{ -\item{matrix}{A dgCMatrix} +\item{x}{A matrix.} } \description{ Attempts to turn a dgCMatrix into a dense matrix } +\examples{ +data("pancreas_sub") +system.time(mat1 <- as.matrix(slot(pancreas_sub[["RNA"]], "counts"))) +system.time(mat2 <- as_matrix(slot(pancreas_sub[["RNA"]], "counts"))) +identical(mat1, mat2) + +} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 94c018ca..ec7932f3 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -11,7 +11,7 @@ Rcpp::Rostream& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get(); #endif // asMatrix -IntegerMatrix asMatrix(NumericVector rp, NumericVector cp, NumericVector z, int nrows, int ncols); +NumericMatrix asMatrix(NumericVector rp, NumericVector cp, NumericVector z, int nrows, int ncols); RcppExport SEXP _SCP_asMatrix(SEXP rpSEXP, SEXP cpSEXP, SEXP zSEXP, SEXP nrowsSEXP, SEXP ncolsSEXP) { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; diff --git a/src/asMatrix.cpp b/src/asMatrix.cpp index 5f9b2ca0..0897a1ba 100644 --- a/src/asMatrix.cpp +++ b/src/asMatrix.cpp @@ -1,21 +1,20 @@ -#include -using namespace Rcpp; +#include +using namespace Rcpp; - -// [[Rcpp::export]] -IntegerMatrix asMatrix(NumericVector rp, +// [[Rcpp::export]] +NumericMatrix asMatrix(NumericVector rp, NumericVector cp, NumericVector z, int nrows, int ncols){ - int k = z.size() ; + int k = z.size() ; - IntegerMatrix mat(nrows, ncols); + NumericMatrix mat(nrows, ncols); - for (int i = 0; i < k; i++){ - mat(rp[i],cp[i]) = z[i]; - } + for (int i = 0; i < k; i++){ + mat(rp[i],cp[i]) = z[i]; + } - return mat; -} \ No newline at end of file + return mat; +} From f0a0b4e7a9039a2c0534537b92a0e1c34daa1346 Mon Sep 17 00:00:00 2001 From: zhanghao-njmu <542370159@qq.com> Date: Tue, 14 Nov 2023 15:13:20 +0800 Subject: [PATCH 4/4] Increment version number to 0.5.5 --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 55efbad8..d22293d1 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: SCP Type: Package Title: Single Cell Pipeline -Version: 0.5.4 +Version: 0.5.5 Author: Hao Zhang Maintainer: Hao Zhang Description: An end-to-end Single-Cell Pipeline designed to facilitate comprehensive analysis and exploration of single-cell data.