diff --git a/.Rbuildignore b/.Rbuildignore index f47a209..cec6fa5 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -2,3 +2,5 @@ ^scratch$ ^test_obj\.rds$ ^utils$ +^\.github$ +^_pkgdown\.yml$ diff --git a/.github/workflows/r-tests.yml b/.github/workflows/r-tests.yml new file mode 100644 index 0000000..9996136 --- /dev/null +++ b/.github/workflows/r-tests.yml @@ -0,0 +1,21 @@ +name: R tests + +on: + pull_request: + +jobs: + test: + runs-on: ubuntu-latest + + steps: + - uses: actions/checkout@v4 + + - uses: r-lib/actions/setup-r@v2 + + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + extra-packages: any::rcmdcheck + needs: test + + - name: Run rcmdcheck + run: Rscript -e 'rcmdcheck::rcmdcheck(args = c("--no-manual"), error_on = "error")' diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..ecf242d --- /dev/null +++ b/.gitignore @@ -0,0 +1,3 @@ +docs +*.csv +*.json \ No newline at end of file diff --git a/DESCRIPTION b/DESCRIPTION index b96c11c..b2872f7 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: AzimuthAPI Title: Pan-Azimuth Web API Interface -Version: 0.2.0 +Version: 0.9.0 Authors@R: person("Satija", "Lab", email = "satijalabnygc@gmail.com", role = c("aut", "cre")) Description: An R package providing an interface to the Pan-Azimuth Web API for single-cell RNA sequencing analysis. @@ -8,20 +8,26 @@ License: MIT + file LICENSE Encoding: UTF-8 Roxygen: list(markdown = TRUE) RoxygenNote: 7.3.3 +VignetteBuilder: knitr Imports: Seurat, SeuratObject, + argparse, + cli, + curl, dplyr, ggplot2, - tibble, httr, jsonlite, - curl, - RCurl, - reticulate, Matrix, - argparse, - yaml, methods, + tibble, + RCurl, + reticulate, + rlang, stats, utils +Suggests: + knitr, + rmarkdown, + testthat diff --git a/NAMESPACE b/NAMESPACE index 9b46d9b..ead5f23 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -3,15 +3,8 @@ export(ANNotate) export(CloudAzimuth) export(PrepLabel) -export(get_data) -export(if_gpu) -export(listen_to_progress) export(make_QC_heatmap) export(make_azimuth_QC_heatmaps) -export(package_obj) -export(process_azimuth_results) -export(read_obj_min) -export(run_azimuth_api) importFrom(Matrix,Matrix) importFrom(RCurl,url.exists) importFrom(Seurat,"DefaultAssay<-") @@ -35,6 +28,7 @@ importFrom(argparse,ArgumentParser) importFrom(curl,curl_fetch_stream) importFrom(curl,form_file) importFrom(curl,handle_setform) +importFrom(curl,handle_setopt) importFrom(curl,new_handle) importFrom(dplyr,"%>%") importFrom(dplyr,across) @@ -54,6 +48,7 @@ importFrom(httr,GET) importFrom(httr,POST) importFrom(httr,content) importFrom(httr,status_code) +importFrom(httr,timeout) importFrom(httr,upload_file) importFrom(httr,write_disk) importFrom(jsonlite,fromJSON) @@ -61,6 +56,8 @@ importFrom(methods,"slot<-") importFrom(methods,slot) importFrom(methods,slotNames) importFrom(reticulate,r_to_py) +importFrom(rlang,"%||%") +importFrom(rlang,.data) importFrom(stats,cor) importFrom(stats,dist) importFrom(stats,hclust) diff --git a/R/annotation.R b/R/annotation.R index e1a19eb..7d9a04b 100644 --- a/R/annotation.R +++ b/R/annotation.R @@ -1,4 +1,8 @@ -#' Run Pan-Human Azimuth annotation +#' Run Pan-human Azimuth annotation locally +#' +#' This function runs the Pan-human Azimuth model on a Seurat object to annotate cell types, via reticulate and the `panhumanpy` Python package. **We recommend using the `CloudAzimuth` function, which runs cloud-based annotation, can handle large datasets, and performs robust error handling.** +#' +#' This function requires the `panhumanpy` Python package to be installed and accessible via reticulate. #' #' @param query_obj Seurat object to annotate #' @param feature_names_col Column name for feature names @@ -8,8 +12,10 @@ #' @param norm_check_batch_size Batch size to inspect normalization of data #' @param output_mode Output mode for annotated cell metadata #' @param refine_labels Whether to refine labels -#' @param extract_embeddings Whether to azimuth embeddings -#' @param umap_embeddings Whether to include UMAP embeddings +#' @param map_to_cl One or more annotation columns to map to Cell Ontology labels +#' @param include_cl_id Whether to add Cell Ontology IDs to the output metadata +#' @param extract_embeddings Whether to extract Azimuth embeddings +#' @param umap_embeddings Whether to include UMAP of Azimuth embeddings #' @param n_neighbors Number of neighbors for UMAP #' @param n_components Number of components for UMAP #' @param metric Distance metric for UMAP @@ -23,10 +29,13 @@ #' @param cutoff_abs Absolute cutoff for label filtering #' @param cutoff_frac Fractional cutoff for label filtering #' @param model_version Version of the model to use +#' @param assay Assay to use for annotation #' #' @importFrom SeuratObject Idents<- #' @importFrom reticulate r_to_py +#' @importFrom rlang %||% #' +#' @concept annotation #' @return Annotated Seurat object #' #' @export @@ -34,30 +43,36 @@ ANNotate <- function( query_obj, feature_names_col = NULL, annotation_pipeline = 'supervised', - eval_batch_size = 40960, + eval_batch_size = 8192L, normalization_override = FALSE, - norm_check_batch_size = 1000, + norm_check_batch_size = 100L, output_mode = 'minimal', refine_labels = TRUE, + map_to_cl = NULL, + include_cl_id = FALSE, extract_embeddings = TRUE, umap_embeddings = TRUE, - n_neighbors = 30, - n_components = 2, + n_neighbors = 30L, + n_components = 2L, metric = "cosine", min_dist = 0.3, umap_lr = 1.0, - umap_seed = 42, + umap_seed = 42L, spread = 1.0, verbose = TRUE, - model_version = "v0", init = "spectral", + model_version = "v0", process_obj = TRUE, cutoff_abs = 5, - cutoff_frac = 0.001 + cutoff_frac = 0.001, + assay = NULL ) { options(warn = -1) - #source_data_dir <- paste0(python_module_path, source_data_dir) - cat("Running Pan-Human Azimuth:\n") + # python dependencies + annotate <- reticulate::import("panhumanpy.ANNotate") + sp <- reticulate::import("scipy.sparse") + + cat("Running Pan-human Azimuth:\n") cat("\n") # Convert integers @@ -65,13 +80,16 @@ ANNotate <- function( n_neighbors <- as.integer(n_neighbors) n_components <- as.integer(n_components) umap_seed <- as.integer(umap_seed) + + assay <- assay %||% DefaultAssay(query_obj) # Read and process the Seurat object - query <- read_obj_min(query_obj, feature_names_col) - X_query <- sp$csr_matrix(r_to_py(query$X_query)) + query <- read_obj_min(query_obj, feature_names_col, assay_default = assay) + X_query <- sp$csr_matrix(reticulate::r_to_py(query$X_query)) query_features <- query$query_features cells_meta <- query$query_cells_df - + assay_cells <- query$assay_cells + # Run annotation core core_outputs <- annotate$annotate_core( X_query, @@ -83,6 +101,8 @@ ANNotate <- function( norm_check_batch_size, output_mode, refine_labels, + map_to_cl, + include_cl_id, extract_embeddings, umap_embeddings, n_neighbors, @@ -93,27 +113,56 @@ ANNotate <- function( umap_seed, spread, verbose, - model_version, - init + init, + model_version ) - embeddings_mode <- core_outputs[[3]] - embeddings_dict <- core_outputs[[4]] - query_cells_df <- core_outputs[[10]] - if_umap_embeddings <- core_outputs[[11]] - umap_embeddings_dict <- core_outputs[[12]] + embeddings_dict <- core_outputs$embeddings_dict + umap_embeddings_dict <- core_outputs$umap_dict + cells_meta <- core_outputs$cells_meta - annotated_obj <- package_obj(embeddings_mode, embeddings_dict, if_umap_embeddings, - umap_embeddings_dict, query_cells_df, query_obj) + annotated_obj = package_obj( + extract_embeddings, + embeddings_dict, + umap_embeddings, + umap_embeddings_dict, + cells_meta, + assay_cells, + query_obj + ) - if (process_obj) { - annotated_obj <- PrepLabel(annotated_obj, - label_id = 'final_level_label', - cutoff = min(cutoff_abs, cutoff_frac * ncol(annotated_obj)), - cutid = 'Other', - newid = 'azimuth_label') + if (process_obj){ + annotated_obj <- PrepLabel(annotated_obj, + label_id = 'final_level_labels', + cutoff = min(cutoff_abs, cutoff_frac * ncol(annotated_obj)), + cutid = 'Other', + newid = 'azimuth_label') Idents(annotated_obj) <- 'azimuth_label' } return(annotated_obj) } + +#' Prepare labels for annotation +#' +#' @param object Seurat object +#' @param label_id Column name for labels +#' @param newid New column name for processed labels +#' @param cutid Label for rejected cells +#' @param cutoff Minimum count threshold +#' @return Updated Seurat object +#' @concept annotation +#' @export +PrepLabel <- function( + object, + label_id = 'final_level_label', + newid = 'PrepLabel', + cutid = 'Other', + cutoff=10 + ) { + rejected_names <- names(which(table(object@meta.data[,label_id])= 300) { + stop("Failed to download output file from the server.") + } } -#' Run the PanAzimuth API analysis +#' Check local vs latest API version / ensure connection can be established to the server #' -#' @param input_data Input data for analysis -#' @param api_endpoint API endpoint URL -#' @return Analysis results -#' @export -run_azimuth_api <- function(input_data, api_endpoint) { - # Implementation of the API analysis - # This is a placeholder - implement the actual API interaction logic - stop("Not implemented yet") -} - -#' Process results from PanAzimuth API +#' @param api_base_url Base URL for the API +#' +#' @importFrom httr GET content status_code timeout +#' @importFrom utils packageVersion #' -#' @param results Raw results from the API -#' @return Processed results in Seurat object format -#' @export -process_azimuth_results <- function(results) { - # Implementation of results processing - # This is a placeholder - implement the actual results processing logic - stop("Not implemented yet") -} \ No newline at end of file +#' @return NULL +#' @keywords internal +#' @noRd +check_api_version <- function(api_base_url) { + version_url <- paste0(api_base_url, "/version") + response <- tryCatch({ + GET(version_url, timeout(5)) + }, error = function(e) { + if (grepl("Could not connect to server", e$message, fixed = TRUE)) { + stop(simpleError( + "Connection refused: server not running or port closed.\nPlease report at https://github.com/satijalab/AzimuthAPI/issues.", + call = conditionCall(e) + )) + } + stop(simpleError(paste0("Error connecting to the server: ", conditionMessage(e)), call = conditionCall(e))) + }) + + if (status_code(response) != 200) { + stop(simpleError(paste0("Failed to retrieve version information from the server (HTTP ", status_code(response), ")."), call = NULL)) + } + + latest_version <- content(response)$version + current_version <- as.character(packageVersion("AzimuthAPI")) + return(utils::compareVersion(current_version, latest_version) < 0) +} \ No newline at end of file diff --git a/R/argument_parser.R b/R/argument_parser.R index e8d3b6b..ce2926f 100644 --- a/R/argument_parser.R +++ b/R/argument_parser.R @@ -6,13 +6,14 @@ #' \dontrun{ #' args <- parse_annotate_args() #' } +#' @noRd parse_annotate_args <- function() { # Check if argparse is available if (!requireNamespace("argparse", quietly = TRUE)) { stop("The 'argparse' package must be installed to parse command line arguments.") } - message("Parsing command line arguments...") + cli::cli_alert_info("Parsing command line arguments...") # Create argument parser parser <- ArgumentParser(description = "Command line interface for ANNotate function") @@ -41,22 +42,22 @@ parse_annotate_args <- function() { parser$add_argument( "-ebs", "--eval_batch_size", - default = 40960, - help = "Batch size for evaluation (default: 40960)", + default = 8192, + help = "Batch size for evaluation (default: 8192)", type = "integer" ) parser$add_argument( "-no", "--normalization_override", - action = "store_true", + action = "store_false", default = FALSE, help = "Override normalization (default: FALSE)" ) parser$add_argument( "-ncb", "--norm_check_batch_size", - default = 1000, - help = "Batch size for normalization check (default: 1000)", + default = 100, + help = "Batch size for normalization check (default: 100)", type = "integer" ) @@ -73,6 +74,20 @@ parse_annotate_args <- function() { default = TRUE, help = "Refine labels (default: TRUE)" ) + + parser$add_argument( + "-mtc", "--map_to_cl", + default = NULL, + help = "Colname(s) to map to Cell Ontology labels (default: NULL)", + type = "character" + ) + + parser$add_argument( + "-icl", "--include_cl_id", + action = "store_false", + default = FALSE, + help = "Include Cell Ontology IDs (default: FALSE)" + ) parser$add_argument( "-ee", "--extract_embeddings", @@ -196,8 +211,9 @@ parse_annotate_args <- function() { #' args <- parse_annotate_args() #' formatted_args <- format_annotate_args(args) #' } +#' @noRd format_annotate_args <- function(args) { - message("Formatting arguments for ANNotate function...") + cli::cli_alert_info("Formatting arguments for ANNotate function...") # Convert arguments to properly named list for ANNotate formatted_args <- list( @@ -212,6 +228,8 @@ format_annotate_args <- function(args) { norm_check_batch_size = args$norm_check_batch_size, output_mode = args$output_mode, refine_labels = args$refine_labels, + map_to_cl = args$map_to_cl, + include_cl_id = args$include_cl_id, extract_embeddings = args$extract_embeddings, umap_embeddings = args$umap_embeddings, n_neighbors = args$n_neighbors, @@ -222,8 +240,8 @@ format_annotate_args <- function(args) { umap_seed = args$umap_seed, spread = args$spread, verbose = args$verbose, - model_version = args$model_version, init = args$init, + model_version = args$model_version, process_obj = args$process_obj, cutoff_abs = args$cutoff_abs, cutoff_frac = args$cutoff_frac @@ -231,15 +249,15 @@ format_annotate_args <- function(args) { # Inform about parsed arguments if (args$verbose) { - message("Arguments for ANNotate function:") + cli::cli_alert_info("Arguments for ANNotate function:") for (name in names(formatted_args)) { value <- formatted_args[[name]] if (is.null(value)) { - message(sprintf(" %s: NULL", name)) + cli::cli_text(" {name}: NULL") } else if (is.logical(value) || is.numeric(value)) { - message(sprintf(" %s: %s", name, value)) + cli::cli_text(" {name}: {value}") } else { - message(sprintf(" %s: '%s'", name, value)) + cli::cli_text(" {name}: '{value}'") } } } diff --git a/R/cloud.R b/R/cloud.R index 146db9e..c411c81 100644 --- a/R/cloud.R +++ b/R/cloud.R @@ -1,25 +1,34 @@ -#' Run Pan-Human Azimuth annotation on the cloud +#' Run Pan-human Azimuth annotation on the cloud #' #' @param object Seurat object to annotate #' @param assay Name of the assay to use (default: 'RNA') #' @param ip Hostname or IP address of the cloud server (default: 'azimuthapi.satijalab.org') #' @param port Port number for the API (default: 5000) -#' @param ... Additional arguments for the API to pass to the model +#' @param ... Additional arguments for the API to pass to the model (see ANNotate function for details) #' @return Annotated Seurat object #' @importFrom httr POST GET upload_file content status_code #' @importFrom RCurl url.exists #' @importFrom SeuratObject LayerData Idents IsMatrixEmpty CreateAssay5Object CreateSeuratObject Cells Idents<- +#' @concept annotation #' @export CloudAzimuth <- function(object = object, assay = 'RNA', ip = 'azimuthapi.satijalab.org', port = 5000, ...) { - message("Running Pan-Human Azimuth on the cloud!") - + + cli::cli_h1("Running Pan-human Azimuth on the cloud") + + api_base_url <- paste0('http://', ip, ":", port) + update <- check_api_version(api_base_url) + + if (isTRUE(update)) { + cli::cli_alert_warning("A new version of the AzimuthAPI package is available. Please update to the latest version for the best experience.") + } + layer_name <- 'data' tryCatch({ data <- LayerData(object, assay = assay, layer = layer_name) }, warning = function(w) { - message <- conditionMessage(w) - if (grepl(paste0("Layer.*", layer_name, ".*is empty"), message)) { + warn_msg <- conditionMessage(w) + if (grepl(paste0("Layer.*", layer_name, ".*is empty"), warn_msg)) { stop(simpleError( "Please run NormalizeData on the data before running Azimuth", call = conditionCall(w) @@ -40,17 +49,15 @@ CloudAzimuth <- function(object = object, assay = 'RNA', ip = 'azimuthapi.satija features <- readLines(feature_file) data <- data[intersect(features, rownames(data)), ] } - - message("Uploading dataset") + suppressWarnings(srt <- CreateSeuratObject(CreateAssay5Object(data = data))) - + tmpname <- tempfile() tmp_input <- paste0(tmpname, ".rds") tmp_output <- paste0(tmpname, "_ANN.rds") saveRDS(object = srt, file = tmp_input) - api_base_url <- paste0('http://', ip, ":", port) - check_api_version(api_base_url) + cli::cli_alert_info("Uploading dataset...") process_rds_file(api_base_url, tmp_input, ...) srt <- readRDS(file = tmp_output) @@ -76,6 +83,10 @@ CloudAzimuth <- function(object = object, assay = 'RNA', ip = 'azimuthapi.satija if (!i %in% colnames(object@meta.data)) { object@meta.data[, i] <- NA } + # Convert to character if it's a factor to avoid level mismatch issues + if (is.factor(object@meta.data[[i]])) { + object@meta.data[[i]] <- as.character(object@meta.data[[i]]) + } # Update only the processed cells object@meta.data[assay_cells, i] <- srt@meta.data[assay_cells, i] } @@ -84,39 +95,6 @@ CloudAzimuth <- function(object = object, assay = 'RNA', ip = 'azimuthapi.satija return(object) } -#' Check local vs latest API version / ensure connection can be established to the server -#' -#' @param api_base_url Base URL for the API - -#' @importFrom httr GET content status_code -#' @importFrom utils packageVersion -#' -#' @return NULL -check_api_version <- function(api_base_url) { - version_url <- paste0(api_base_url, "/version") - tryCatch({ - response <- GET(version_url) - }, error = function(e) { - if (grepl("Could not connect to server", e$message, fixed = TRUE)) { - stop(simpleError( - "Connection refused: server not running or port closed.\nPlease report at https://github.com/satijalab/AzimuthAPI/issues.", - call = conditionCall(e) - )) - } - }) - - if (status_code(response) == 200) { - version_info <- content(response) - version <- version_info$version - if (packageVersion("AzimuthAPI") < version) { - warning("A new version of the AzimuthAPI package is available: ", - version, "\n", "Please update to the latest version for the best experience.\n") - } - } else { - warning("Failed to retrieve API version information.\n") - } -} - #' Process RDS file through the cloud API #' #' @param api_base_url Base URL for the API @@ -124,20 +102,35 @@ check_api_version <- function(api_base_url) { #' @param ... Additional arguments passed to the API #' @return NULL #' @importFrom httr POST GET upload_file content status_code +#' @noRd process_rds_file <- function(api_base_url, file_path, ...) { progress_url <- paste0(api_base_url, "/process_rds") - cat("Uploading file and listening for updates...\n") - success <- listen_to_progress(progress_url, file_path, ...) + cli::cli_alert_info("Uploading file and listening for updates...") + result <- listen_to_progress(progress_url, file_path, ...) - if (isFALSE(success)) { + if (isFALSE(result$success)) { stop("Processing failed on the server. Please check the error messages above.") } - + output_file_name <- gsub("\\.rds$", "_ANN.rds", basename(file_path)) - download_url <- paste0(api_base_url, "/download_output?output_file=/tmp/", - output_file_name) + if (!is.null(result$output_file)) { + output_file_name <- basename(result$output_file) + } + download_url <- result$download_url + if (is.null(download_url) && !is.null(result$output_file)) { + download_url <- paste0(api_base_url, "/download_output?output_file=", + utils::URLencode(result$output_file, reserved = TRUE) + ) + } + if (is.null(download_url)) { + stop("Server did not return output file metadata for download.") + } + if (!grepl("^https?://", download_url)) { + download_url <- paste0(api_base_url, download_url) + } save_path <- file.path(tempdir(), output_file_name) - cat("Downloading the output file...\n") + cli::cli_alert_info("Downloading the output file...") download_output(download_url, save_path) + cli::cli_alert_success("Annotation complete. Output saved to: {save_path}") } diff --git a/R/gpu_utils.R b/R/gpu_utils.R index 800a228..24c705a 100644 --- a/R/gpu_utils.R +++ b/R/gpu_utils.R @@ -1,7 +1,8 @@ #' Check for NVIDIA GPU availability #' #' @return Logical indicating whether an NVIDIA GPU is available -#' @export +#' @keywords internal +#' @noRd if_gpu <- function() { if (.Platform$OS.type == "windows") { # Windows-specific check diff --git a/R/make_qc_heatmaps.R b/R/make_qc_heatmaps.R index c021635..3a1c2df 100644 --- a/R/make_qc_heatmaps.R +++ b/R/make_qc_heatmaps.R @@ -18,7 +18,7 @@ #' #' @importFrom Seurat FindAllMarkers ScaleData DoHeatmap SetIdent NoLegend #' @importFrom SeuratObject Idents Idents<- LayerData Cells -#' @importFrom dplyr group_by filter slice_head ungroup left_join summarise across %>% mutate where arrange +#' @importFrom dplyr group_by filter slice_head ungroup left_join summarise across mutate where arrange #' @importFrom stats dist cor hclust #' @importFrom tibble rownames_to_column column_to_rownames #' @importFrom ggplot2 theme element_text ggsave @@ -51,8 +51,8 @@ make_QC_heatmap <- function( } if (!is.null(n_downsample)) seurat_obj <- subset(seurat_obj, downsample=n_downsample) mark_all <- FindAllMarkers(seurat_obj,only.pos = TRUE,min.pct = min.pct) - mark_all %>% group_by(cluster) %>% - dplyr::filter(avg_log2FC > logfc_cutoff) %>% + mark_all %>% group_by(.data$cluster) %>% + dplyr::filter(.data$avg_log2FC > logfc_cutoff) %>% slice_head(n = n_markers) %>% ungroup() -> top_markers @@ -70,7 +70,7 @@ make_QC_heatmap <- function( as.data.frame() %>% rownames_to_column("Cell") %>% left_join(data.frame(Cell = Cells(seurat_obj), CellType = Idents(seurat_obj)), by = "Cell") %>% - group_by(CellType) %>% + group_by(.data$CellType) %>% summarise(across(where(is.numeric), mean, na.rm = TRUE)) %>% # Use `where(is.numeric)` to avoid non-numeric columns column_to_rownames("CellType") @@ -86,8 +86,8 @@ make_QC_heatmap <- function( # Step 5: Reorder the top_markers dataframe based on ordered_cell_types top_markers <- top_markers %>% - mutate(cluster = factor(cluster, levels = ordered_cell_types)) %>% - arrange(cluster) + mutate(cluster = factor(.data$cluster, levels = ordered_cell_types)) %>% + arrange(.data$cluster) } cells.plot <- (names(which(!is.na(Idents(seurat_obj))))) @@ -130,7 +130,7 @@ make_azimuth_QC_heatmaps <- function( # filter non-concordant # SKYLAR: Proposed fix due to ANNotate metadata column name change - metadata <- subset(metadata, full_consistent_hierarchy==TRUE) + metadata <- metadata[metadata$full_consistent_hierarchy==TRUE, ] level2_name = 'level_two_labels' metadata[, level2_name] <- sapply(strsplit(metadata[, full_name], "\\|"), function(x) ifelse(length(x) > 1, x[2], "")) @@ -190,7 +190,7 @@ make_azimuth_QC_heatmaps <- function( tryCatch({ plot_list[[level1]] <- make_QC_heatmap(lobj, min.size = min.final.group, identity = as.character(level1), ...) }, error = function(e) { - message(paste("Error in processing", level1, ":", e$message)) + cli::cli_warn("Error in processing {level1}: {e$message}") }) } return(plot_list) diff --git a/R/read_seurat_object.R b/R/read_seurat_object.R deleted file mode 100644 index add3fb8..0000000 --- a/R/read_seurat_object.R +++ /dev/null @@ -1,58 +0,0 @@ -#' Read a Seurat object from an RDS file -#' -#' Reads a Seurat object from an RDS file, ensuring -#' compatibility with Seurat versions >= 4.4.0 -#' -#' @param filepath Character string specifying the path to the RDS file -#' @param assay_name Character string specifying the assay to use (default: "RNA") -#' @importFrom utils packageVersion -#' @importFrom Seurat DefaultAssay<- -#' -#' @return A valid Seurat object -#' -#' @examples -#' \dontrun{ -#' seu_obj <- read_seurat_object("path/to/seurat_object.rds") -#' } -read_seurat_object <- function(filepath, assay_name = "RNA") { - # Check if file exists - if (!file.exists(filepath)) { - stop(paste("File does not exist at the specified path:", filepath)) - } - - # Check file extension - if (!grepl("\\.rds$", filepath, ignore.case = TRUE)) { - warning("File does not have .rds extension. Attempting to read anyway.") - } - - # Load the RDS file - tryCatch({ - seurat_obj <- readRDS(filepath) - }, error = function(e) { - stop(paste("Error reading RDS file:", e$message)) - }) - - # Check if it's a Seurat object - if (!inherits(seurat_obj, "Seurat")) { - stop("The file does not contain a Seurat object.") - } - - # Check Seurat version - seurat_version <- packageVersion("Seurat") - if (seurat_version < "4.4.0") { - warning(paste("Current Seurat version:", seurat_version, - "is below 4.4.0. Some functionality may not work as expected.")) - } - - # Check if specified assay exists - if (!assay_name %in% names(seurat_obj@assays)) { - available_assays <- paste(names(seurat_obj@assays), collapse = ", ") - stop(paste("Assay", assay_name, "not found in the Seurat object.", - "Available assays:", available_assays)) - } - - # Set default assay - DefaultAssay(seurat_obj) <- assay_name - - return(seurat_obj) -} \ No newline at end of file diff --git a/R/save_seurat_object.R b/R/save_seurat_object.R deleted file mode 100644 index 7b09bba..0000000 --- a/R/save_seurat_object.R +++ /dev/null @@ -1,78 +0,0 @@ -#' Save a Seurat object to an RDS file -#' -#' Saves a Seurat object to an RDS file after applying "keto diet" to reduce size -#' -#' @param seurat_obj A Seurat object to save -#' @param filepath Character string specifying the output path -#' @importFrom utils packageVersion -#' @importFrom methods slotNames slot slot<- -#' -#' @return Invisible filepath where the object was saved -#' -save_seurat_object <- function(seurat_obj, filepath) { - # Check if object is a Seurat object - if (!inherits(seurat_obj, "Seurat")) { - stop("The object is not a Seurat object.") - } - - # Construct output filepath - dir_path <- dirname(filepath) - file_name <- basename(filepath) - file_name_ann <- sub("\\.rds$", "_ANN.rds", file_name) - if (file_name == file_name_ann) { - # If no .rds extension was found, add it - file_name_ann <- paste0(file_name, "_ANN.rds") - } - output_path <- file.path(dir_path, file_name_ann) - - # This is a joke about massively dieting the object - # so that we don't redownload expression or embedding data - # in this API case we know the object has an RNA assay with counts and data layers - keto_object=TRUE - if (keto_object) { - if ("RNA" %in% names(seurat_obj@assays)) { - # Handle Seurat v5+ and v4 differently - if (packageVersion("Seurat") >= "5.0.0") { - # For Seurat v5+ - if ("data" %in% names(seurat_obj[["RNA"]]@layers)) { - seurat_obj[["RNA"]]$data <- Matrix::Matrix(0, - nrow = nrow(seurat_obj[["RNA"]]$data), - ncol = ncol(seurat_obj[["RNA"]]$data), - sparse = TRUE) - } - if ("counts" %in% names(seurat_obj[["RNA"]]@layers)) { - seurat_obj[["RNA"]]$counts <- Matrix::Matrix(0, - nrow = nrow(seurat_obj[["RNA"]]$counts), - ncol = ncol(seurat_obj[["RNA"]]$counts), - sparse = TRUE) - } - } else { - # For Seurat v4 - if ("data" %in% slotNames(seurat_obj[["RNA"]])) { - slot(seurat_obj[["RNA"]], "data") <- Matrix::Matrix(0, - nrow = nrow(slot(seurat_obj[["RNA"]], "data")), - ncol = ncol(slot(seurat_obj[["RNA"]], "data")), - sparse = TRUE) - } - if ("counts" %in% slotNames(seurat_obj[["RNA"]])) { - slot(seurat_obj[["RNA"]], "counts") <- Matrix::Matrix(0, - nrow = nrow(slot(seurat_obj[["RNA"]], "counts")), - ncol = ncol(slot(seurat_obj[["RNA"]], "counts")), - sparse = TRUE) - } - } - } - - # Remove ANNshallow_embeddings - if ("azimuth_embed" %in% names(seurat_obj@reductions)) { - seurat_obj[["azimuth_embed"]] <- NULL - } - } - - # Save the object - message("Saving annotated object to ", output_path) - saveRDS(seurat_obj, file = output_path) - - # Return the path invisibly - invisible(output_path) -} \ No newline at end of file diff --git a/R/seurat_utils.R b/R/seurat_utils.R index bf2fafc..01deced 100644 --- a/R/seurat_utils.R +++ b/R/seurat_utils.R @@ -1,3 +1,151 @@ +#' @importFrom rlang %||% .data +#' @importFrom dplyr %>% +NULL + + +#' Read a Seurat object from an RDS file +#' +#' Reads a Seurat object from an RDS file, ensuring +#' compatibility with Seurat versions >= 4.4.0 +#' +#' @param filepath Character string specifying the path to the RDS file +#' @param assay_name Character string specifying the assay to use (default: "RNA") +#' @importFrom utils packageVersion +#' @importFrom Seurat DefaultAssay<- +#' +#' @return A valid Seurat object +#' +#' @keywords internal +#' @noRd +#' +#' @examples +#' \dontrun{ +#' seu_obj <- read_seurat_object("path/to/seurat_object.rds") +#' } +read_seurat_object <- function(filepath, assay_name = "RNA") { + # Check if file exists + if (!file.exists(filepath)) { + stop(paste("File does not exist at the specified path:", filepath)) + } + + # Check file extension + if (!grepl("\\.rds$", filepath, ignore.case = TRUE)) { + warning("File does not have .rds extension. Attempting to read anyway.") + } + + # Load the RDS file + tryCatch({ + seurat_obj <- readRDS(filepath) + }, error = function(e) { + stop(paste("Error reading RDS file:", e$message)) + }) + + # Check if it's a Seurat object + if (!inherits(seurat_obj, "Seurat")) { + stop("The file does not contain a Seurat object.") + } + + # Check Seurat version + seurat_version <- packageVersion("Seurat") + if (seurat_version < "4.4.0") { + warning(paste("Current Seurat version:", seurat_version, + "is below 4.4.0. Some functionality may not work as expected.")) + } + + # Check if specified assay exists + if (!assay_name %in% names(seurat_obj@assays)) { + available_assays <- paste(names(seurat_obj@assays), collapse = ", ") + stop(paste("Assay", assay_name, "not found in the Seurat object.", + "Available assays:", available_assays)) + } + + # Set default assay + DefaultAssay(seurat_obj) <- assay_name + + return(seurat_obj) +} + +#' Save a Seurat object to an RDS file +#' +#' Saves a Seurat object to an RDS file after applying "keto diet" to reduce size +#' +#' @param seurat_obj A Seurat object to save +#' @param filepath Character string specifying the output path +#' @importFrom utils packageVersion +#' @importFrom methods slotNames slot slot<- +#' +#' @return Invisible filepath where the object was saved +#' @keywords internal +#' @noRd +#' +save_seurat_object <- function(seurat_obj, filepath) { + # Check if object is a Seurat object + if (!inherits(seurat_obj, "Seurat")) { + stop("The object is not a Seurat object.") + } + + # Construct output filepath + dir_path <- dirname(filepath) + file_name <- basename(filepath) + file_name_ann <- sub("\\.rds$", "_ANN.rds", file_name) + if (file_name == file_name_ann) { + # If no .rds extension was found, add it + file_name_ann <- paste0(file_name, "_ANN.rds") + } + output_path <- file.path(dir_path, file_name_ann) + + # This is a joke about massively dieting the object + # so that we don't redownload expression or embedding data + # in this API case we know the object has an RNA assay with counts and data layers + keto_object=TRUE + if (keto_object) { + if ("RNA" %in% names(seurat_obj@assays)) { + # Handle Seurat v5+ and v4 differently + if (packageVersion("Seurat") >= "5.0.0") { + # For Seurat v5+ + if ("data" %in% names(seurat_obj[["RNA"]]@layers)) { + seurat_obj[["RNA"]]$data <- Matrix::Matrix(0, + nrow = nrow(seurat_obj[["RNA"]]$data), + ncol = ncol(seurat_obj[["RNA"]]$data), + sparse = TRUE) + } + if ("counts" %in% names(seurat_obj[["RNA"]]@layers)) { + seurat_obj[["RNA"]]$counts <- Matrix::Matrix(0, + nrow = nrow(seurat_obj[["RNA"]]$counts), + ncol = ncol(seurat_obj[["RNA"]]$counts), + sparse = TRUE) + } + } else { + # For Seurat v4 + if ("data" %in% slotNames(seurat_obj[["RNA"]])) { + slot(seurat_obj[["RNA"]], "data") <- Matrix::Matrix(0, + nrow = nrow(slot(seurat_obj[["RNA"]], "data")), + ncol = ncol(slot(seurat_obj[["RNA"]], "data")), + sparse = TRUE) + } + if ("counts" %in% slotNames(seurat_obj[["RNA"]])) { + slot(seurat_obj[["RNA"]], "counts") <- Matrix::Matrix(0, + nrow = nrow(slot(seurat_obj[["RNA"]], "counts")), + ncol = ncol(slot(seurat_obj[["RNA"]], "counts")), + sparse = TRUE) + } + } + } + + # Remove ANNshallow_embeddings + if ("azimuth_embed" %in% names(seurat_obj@reductions)) { + seurat_obj[["azimuth_embed"]] <- NULL + } + } + + # Save the object + message("Saving annotated object to ", output_path) + saveRDS(seurat_obj, file = output_path) + + # Return the path invisibly + invisible(output_path) +} + #' Get data from a Seurat object layer #' #' @param object Seurat object @@ -7,7 +155,8 @@ #' @importFrom utils packageVersion #' @importFrom SeuratObject LayerData #' @return Layer data -#' @export +#' @keywords internal +#' @noRd get_data <- function(object, assay, layer = "data") { if (packageVersion("Seurat") >= "5.0.0") { return(LayerData(object[[assay]], layer = layer)) @@ -17,8 +166,6 @@ get_data <- function(object, assay, layer = "data") { } } - - #' Read and process a Seurat object #' #' @param query_obj Seurat object @@ -30,7 +177,8 @@ get_data <- function(object, assay, layer = "data") { #' @importFrom utils packageVersion #' @importFrom methods slotNames #' @return List containing processed data -#' @export +#' @keywords internal +#' @noRd read_obj_min <- function(query_obj, feature_names_col, assay_default='RNA') { if (!(assay_default %in% names(query_obj@assays))){ @@ -51,12 +199,12 @@ read_obj_min <- function(query_obj, feature_names_col, assay_default='RNA') { query_obj <- NormalizeData(query_obj[[assay_default]]) normalized_data <- get_data(query_obj, assay = assay_default) } - - X_query <- t(normalized_data) + X_query <- Matrix::t(normalized_data) X_query <- Matrix(X_query, sparse = TRUE) + assay_cells <- colnames(normalized_data) cell_metadata <- query_obj@meta.data - query_cells_df <- as.data.frame(cell_metadata) + query_cells_df <- as.data.frame(cell_metadata[assay_cells, , drop = FALSE]) if (!is.null(feature_names_col)) { feature_metacols <- colnames(query_obj[[assay_default]][[]]) @@ -77,7 +225,8 @@ read_obj_min <- function(query_obj, feature_names_col, assay_default='RNA') { return(list( X_query = X_query, query_features = query_features, - query_cells_df = query_cells_df + query_cells_df = query_cells_df, + assay_cells = assay_cells )) } @@ -90,33 +239,47 @@ read_obj_min <- function(query_obj, feature_names_col, assay_default='RNA') { #' @param umap_embeddings Whether UMAP embeddings are computed #' @param umap_embeddings_dict Dictionary of UMAP embeddings #' @param query_cells_df Cell metadata +#' @param query_cells Cells present in the processed assay #' @param query_obj Seurat object #' @importFrom SeuratObject CreateDimReducObject #' @importFrom Seurat DefaultAssay #' @return Updated Seurat object -#' @export +#' @keywords internal +#' @noRd package_obj <- function( extract_embeddings, embeddings_dict, umap_embeddings, umap_embeddings_dict, query_cells_df, + query_cells, query_obj ) { + if (nrow(query_cells_df) != length(query_cells)) { + stop("Dimension mismatch: query cells does not match the processed assay cells.") + } + + if (!identical(rownames(query_cells_df), query_cells)) { + if (setequal(rownames(query_cells_df), query_cells)) { + query_cells_df <- query_cells_df[query_cells, , drop = FALSE] + } else { + rownames(query_cells_df) <- query_cells + } + } if (extract_embeddings) { for (em_name in names(embeddings_dict)) { em_matrix <- as.matrix(embeddings_dict[[em_name]]) - if (nrow(em_matrix) != length(Cells(query_obj))) { + if (nrow(em_matrix) != length(query_cells)) { stop(paste( "Dimension mismatch:", em_name, " does not have as many ", - "cells as the query obj." + "cells as the processed assay." )) } - rownames(em_matrix) <- Cells(query_obj) + rownames(em_matrix) <- query_cells dimreduc_obj <- CreateDimReducObject( embeddings = em_matrix, @@ -133,14 +296,14 @@ package_obj <- function( em_matrix <- as.matrix(umap_embeddings_dict[[em_name]]) - if (nrow(em_matrix) != length(Cells(query_obj))) { + if (nrow(em_matrix) != length(query_cells)) { stop(paste( "Dimension mismatch:", em_name, " does not ", - "have as many cells as the query obj." + "have as many cells as the processed assay." )) } - rownames(em_matrix) <- Cells(query_obj) + rownames(em_matrix) <- query_cells dimreduc_obj <- CreateDimReducObject( embeddings = em_matrix, @@ -152,31 +315,12 @@ package_obj <- function( } } - query_obj@meta.data <- as.data.frame(query_cells_df) + for (md_col in colnames(query_cells_df)) { + if (!md_col %in% colnames(query_obj@meta.data)) { + query_obj@meta.data[, md_col] <- NA + } + query_obj@meta.data[query_cells, md_col] <- query_cells_df[query_cells, md_col] + } return(query_obj) -} - - -#' Prepare labels for annotation -#' -#' @param object Seurat object -#' @param label_id Column name for labels -#' @param newid New column name for processed labels -#' @param cutid Label for rejected cells -#' @param cutoff Minimum count threshold -#' @return Updated Seurat object -#' @export -PrepLabel <- function( - object, - label_id = 'final_level_label', - newid = 'PrepLabel', - cutid = 'Other', - cutoff=10 - ) { - rejected_names <- names(which(table(object@meta.data[,label_id]) [!TIP] +> For more details, including an introductory vignette & function reference, visit https://satijalab.org/pan_human_azimuth. -You can install the package using devtools: +Two options for annotation are available: + +- `CloudAzimuth`: computation occurs on the cloud; requires no additional setup. +- `ANNotate`: computation occurs entirely locally; requires a working Python installation with [`panhumanpy`](https://pypi.org/project/panhumanpy/) installed. + +## Installation ```r -if (!requireNamespace("devtools", quietly = TRUE)) { - install.packages("devtools") -} -devtools::install_github("satijalab/AzimuthAPI") +# Install devtools if not already installed +install.packages("remotes") + +# Install AzimuthAPI from GitHub +remotes::install_github("satijalab/AzimuthAPI") ``` + diff --git a/_pkgdown.yml b/_pkgdown.yml new file mode 100644 index 0000000..21c1f95 --- /dev/null +++ b/_pkgdown.yml @@ -0,0 +1,3 @@ +template: + bootstrap: 5 + light-switch: true diff --git a/docs/articles/quickstart_vignette.html b/docs/articles/quickstart_vignette.html deleted file mode 100644 index 8b8c285..0000000 --- a/docs/articles/quickstart_vignette.html +++ /dev/null @@ -1,249 +0,0 @@ - - - - - - - -PanAzimuth Tutorial • PanAzimuth - - - - - - - - - - - - -
-
- - - - -
-
-
-

PanAzimuth Tutorial

- -

Compiled: March 26, 2025

- - - - -
- - - -
-

Setup the Environment -

-

PanAzimuthWebAPI can be installed from Github using the following commands:

-
-if (!requireNamespace("devtools", quietly = TRUE)) {
-    install.packages("devtools")
-}
-devtools::install_github("rsatija/PanAzimuthWebAPI")
-
-
-

Run ANNotate -

-

For this vignette, we will demonstrate annotating a dataset of human bone marrow mononuclear (BMNC) cells from eight individual donors produced by the Human Cell Atlas. This dataset is available through SeuratData.

- -
-bmcite <- InstallData("bmcite")
-
-bmcite <- LoadData("bmcite")
-

The CloudANNotate function runs Pan-human Azimuth cell type predictions on a Seurat object via a cloud-based API deployed on AWS, and returns results stored in the object’s cell-level metadata.

-
-bmcite <- CloudANNotate(bmcite, ip = "34.222.135.233")
-
-bmcite
-
## An object of class Seurat 
-## 17034 features across 30672 samples within 2 assays 
-## Active assay: RNA (17009 features, 2000 variable features)
-##  2 layers present: counts, data
-##  1 other assay present: ADT
-##  2 dimensional reductions calculated: spca, azimuth_umap
-
-
-

Visualize predictions -

-

You can visualize the predicted cell type labels on azimuthumap, the UMAP reduction of the embedding that the PanAzimuth model learns on the data based on cell type representation.

-

For example, below are the original cell type annotations visualized on azimuthumap:

-
-p0 <- DimPlot(bmcite, group.by = "celltype.l1", label.size = 3, repel = T, label = T) + NoLegend()
-p1 <- DimPlot(bmcite, group.by = "celltype.l2", label.size = 3, repel = T, label = T) + NoLegend()
-p0 | p1
-

-

Pan-human Azimuth returns predictions in multiple formats, as well as softmax probability scores to estimate model confidence:

-
    -
  • -full_hierarchical_labels: The full predicted label for each cell, where ‘|’ delimits each hierarchical level
    • -
  • -level_zero_labels: The most broad level of the full hierarchical label
    • -
  • -final_level_labels: The most granular level of the full hierarchical label
    • -
  • -final_level_softmax_prob: The model’s predicted probability for each cell’s assigned final_level_label, ranging from 0 to 1
    • -
-

In addition, we include the full_consistent_hierarchy flag to indicate whether if the full hierarchical label is valid according to our internal cell type hierarchy.

-
-# All cell-level metadata that is returned
-colnames(bmcite[[]])
-
##  [1] "orig.ident"                "nCount_RNA"               
-##  [3] "nFeature_RNA"              "nCount_ADT"               
-##  [5] "nFeature_ADT"              "lane"                     
-##  [7] "donor"                     "celltype.l1"              
-##  [9] "celltype.l2"               "RNA.weight"               
-## [11] "full_hierarchical_labels"  "level_zero_labels"        
-## [13] "final_level_labels"        "final_level_softmax_prob" 
-## [15] "full_consistent_hierarchy" "azimuth_broad"            
-## [17] "azimuth_medium"            "azimuth_fine"             
-## [19] "azimuth_label"
-

The full_hierarchical_label for each cell provides highly detailed information, but may be difficult to interpret visually compared to level-wise labels.

-
-p2 <- DimPlot(bmcite, group.by = "full_hierarchical_labels", label.size = 1.5, repel = T, label = T) +
-    NoLegend()
-p2
-

-
-p3 <- DimPlot(bmcite, group.by = "final_level_labels", label.size = 3, repel = T, label = T) + NoLegend()
-p4 <- FeaturePlot(bmcite, features = "final_level_softmax_prob")
-p3 | p4
-

-

We also postprocess our predictions to provide labels at three consistent levels of granularity for easy handling, marking any cell with an invalid full hierarchical label (based on full_consistent_hierarchy) as False.

-
    -
  • -azimuth_broad: Corresponds to level_zero_labels (i.e. Immune cell)
    • -
  • -azimuth_medium: Medium level of granularity (i.e. T cell)
    • -
  • -azimuth_fine: High level of granularity (i.e. Memory CD4 T cell)
    • -
-

For azimuth_medium and azimuth_fine, labels have been adjusted to return labels at different levels in the hierarchy to be consistent with common understandings of cell type granularity. This means that for a certain number of cells, an azimuth_medium and/or azimuth_fine label could be more granular than its initial prediction. In such cases, these labels are likely to be an accurate description of the cell, as we use internal calculations based on probability values returned by the model, but should be handled with caution for more in-depth or formal analysis.

-
-p5 <- DimPlot(bmcite, group.by = "azimuth_medium", label.size = 3, repel = T, label = T) + NoLegend()
-p5
-

-
-p6 <- DimPlot(bmcite, group.by = "azimuth_fine", label.size = 3, repel = T, label = T) + NoLegend()
-p6
-

-

To remove the number of labels displayed, you can filter labels with less than a certain number of cells using PrepLabel.

-
-bmcite <- PrepLabel(bmcite, "azimuth_fine", "azimuth_fine_filtered", cutoff = 20)
-p7 <- DimPlot(bmcite, group.by = "azimuth_fine_filtered", label.size = 3, repel = T, label = T) +
-    NoLegend()
-p7
-

-
-
-

Visualize differentially expressed features -

-

The make_azimuth_QC_heatmaps function allows you to easily explore the quality of predicted labels by creating expression heatmaps by predicted cell type, with optional parameters for improved visualization:

-
    -
  • -final_name: Name of metadata columns to group cells by (default is azimuth_fine)
    • -
  • -min.final.group: Minimum number of cells under a cell type to be displayed
    • -
  • -max.ids.per.plot: Number of cell type labels displayed per plot
    • -
  • -reorder: Flag to indicate whether to reorder cell types by transcriptional similarity
    • -
  • -cells.order: Cell names to specify order of cells (i.e. by softmax probability)
    • -
  • -save_folder_path: Save plots as PNG files under specified folder path
    • -
-

Plots are saved by azimuth_broad categories by default, with the exception of immune cell types grouped separately by lymphoid or myeloid/erythroid subpopulations.

-
-plots <- make_azimuth_QC_heatmaps(bmcite)
-
-p8 <- plots[["Immune_Lymphoid cell_1"]]
-p8
-

-
-p9 <- plots[["Immune_Myeloid cell_1"]]
-p9
-

-
-
- - - -
- - - -
-

-

Developed by Rahul Satija.

-
- -
-

-

Site built with pkgdown 2.0.7.

-
- -
-
- - - - - - - - diff --git a/docs/articles/quickstart_vignette_files/figure-html/lymphoid.heatmap-1.png b/docs/articles/quickstart_vignette_files/figure-html/lymphoid.heatmap-1.png deleted file mode 100644 index e850ace..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/lymphoid.heatmap-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/myeloid.heatmap-1.png b/docs/articles/quickstart_vignette_files/figure-html/myeloid.heatmap-1.png deleted file mode 100644 index df5ac9e..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/myeloid.heatmap-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/orig.labels-1.png b/docs/articles/quickstart_vignette_files/figure-html/orig.labels-1.png deleted file mode 100644 index 9f3e9b1..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/orig.labels-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.1-1.png b/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.1-1.png deleted file mode 100644 index ef8a0ec..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.1-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.2-1.png b/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.2-1.png deleted file mode 100644 index dda6e91..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/panazimuth.labels.2-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/prepped.labels-1.png b/docs/articles/quickstart_vignette_files/figure-html/prepped.labels-1.png deleted file mode 100644 index b7ec784..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/prepped.labels-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/refined.labels.fine-1.png b/docs/articles/quickstart_vignette_files/figure-html/refined.labels.fine-1.png deleted file mode 100644 index fb81f99..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/refined.labels.fine-1.png and /dev/null differ diff --git a/docs/articles/quickstart_vignette_files/figure-html/refined.labels.medium-1.png b/docs/articles/quickstart_vignette_files/figure-html/refined.labels.medium-1.png deleted file mode 100644 index ec0b8dc..0000000 Binary files a/docs/articles/quickstart_vignette_files/figure-html/refined.labels.medium-1.png and /dev/null differ diff --git a/man/ANNotate.Rd b/man/ANNotate.Rd index 1446651..5c3e758 100644 --- a/man/ANNotate.Rd +++ b/man/ANNotate.Rd @@ -2,32 +2,35 @@ % Please edit documentation in R/annotation.R \name{ANNotate} \alias{ANNotate} -\title{Run Pan-Human Azimuth annotation} +\title{Run Pan-human Azimuth annotation locally} \usage{ ANNotate( query_obj, feature_names_col = NULL, annotation_pipeline = "supervised", - eval_batch_size = 40960, + eval_batch_size = 8192L, normalization_override = FALSE, - norm_check_batch_size = 1000, + norm_check_batch_size = 100L, output_mode = "minimal", refine_labels = TRUE, + map_to_cl = NULL, + include_cl_id = FALSE, extract_embeddings = TRUE, umap_embeddings = TRUE, - n_neighbors = 30, - n_components = 2, + n_neighbors = 30L, + n_components = 2L, metric = "cosine", min_dist = 0.3, umap_lr = 1, - umap_seed = 42, + umap_seed = 42L, spread = 1, verbose = TRUE, - model_version = "v0", init = "spectral", + model_version = "v0", process_obj = TRUE, cutoff_abs = 5, - cutoff_frac = 0.001 + cutoff_frac = 0.001, + assay = NULL ) } \arguments{ @@ -47,9 +50,13 @@ ANNotate( \item{refine_labels}{Whether to refine labels} -\item{extract_embeddings}{Whether to azimuth embeddings} +\item{map_to_cl}{One or more annotation columns to map to Cell Ontology labels} + +\item{include_cl_id}{Whether to add Cell Ontology IDs to the output metadata} -\item{umap_embeddings}{Whether to include UMAP embeddings} +\item{extract_embeddings}{Whether to extract Azimuth embeddings} + +\item{umap_embeddings}{Whether to include UMAP of Azimuth embeddings} \item{n_neighbors}{Number of neighbors for UMAP} @@ -67,19 +74,25 @@ ANNotate( \item{verbose}{Whether to show progress} -\item{model_version}{Version of the model to use} - \item{init}{Initialization method for UMAP} +\item{model_version}{Version of the model to use} + \item{process_obj}{Whether to process the object} \item{cutoff_abs}{Absolute cutoff for label filtering} \item{cutoff_frac}{Fractional cutoff for label filtering} + +\item{assay}{Assay to use for annotation} } \value{ Annotated Seurat object } \description{ -Run Pan-Human Azimuth annotation +This function runs the Pan-human Azimuth model on a Seurat object to annotate cell types, via reticulate and the \code{panhumanpy} Python package. \strong{We recommend using the \code{CloudAzimuth} function, which runs cloud-based annotation, can handle large datasets, and performs robust error handling.} +} +\details{ +This function requires the \code{panhumanpy} Python package to be installed and accessible via reticulate. } +\concept{annotation} diff --git a/man/CloudAzimuth.Rd b/man/CloudAzimuth.Rd index 626cfbe..7bcf3c9 100644 --- a/man/CloudAzimuth.Rd +++ b/man/CloudAzimuth.Rd @@ -2,12 +2,12 @@ % Please edit documentation in R/cloud.R \name{CloudAzimuth} \alias{CloudAzimuth} -\title{Run Pan-Human Azimuth annotation on the cloud} +\title{Run Pan-human Azimuth annotation on the cloud} \usage{ CloudAzimuth( object = object, assay = "RNA", - ip = "127.0.0.1", + ip = "azimuthapi.satijalab.org", port = 5000, ... ) @@ -21,11 +21,12 @@ CloudAzimuth( \item{port}{Port number for the API (default: 5000)} -\item{...}{Additional arguments for the API to pass to the model} +\item{...}{Additional arguments for the API to pass to the model (see ANNotate function for details)} } \value{ Annotated Seurat object } \description{ -Run Pan-Human Azimuth annotation on the cloud +Run Pan-human Azimuth annotation on the cloud } +\concept{annotation} diff --git a/man/PrepLabel.Rd b/man/PrepLabel.Rd index 068de5b..f3d600c 100644 --- a/man/PrepLabel.Rd +++ b/man/PrepLabel.Rd @@ -1,5 +1,5 @@ % Generated by roxygen2: do not edit by hand -% Please edit documentation in R/seurat_utils.R +% Please edit documentation in R/annotation.R \name{PrepLabel} \alias{PrepLabel} \title{Prepare labels for annotation} @@ -29,3 +29,4 @@ Updated Seurat object \description{ Prepare labels for annotation } +\concept{annotation} diff --git a/man/check_api_version.Rd b/man/check_api_version.Rd deleted file mode 100644 index 5b3b867..0000000 --- a/man/check_api_version.Rd +++ /dev/null @@ -1,14 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/cloud.R -\name{check_api_version} -\alias{check_api_version} -\title{Check local vs latest API version / ensure connection can be established to the server} -\usage{ -check_api_version(api_base_url) -} -\arguments{ -\item{api_base_url}{Base URL for the API} -} -\description{ -Check local vs latest API version / ensure connection can be established to the server -} diff --git a/man/download_output.Rd b/man/download_output.Rd deleted file mode 100644 index 1214d5d..0000000 --- a/man/download_output.Rd +++ /dev/null @@ -1,16 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/api_interface.R -\name{download_output} -\alias{download_output} -\title{Download output file from the API} -\usage{ -download_output(url, save_path) -} -\arguments{ -\item{url}{Download URL} - -\item{save_path}{Path to save the downloaded file} -} -\description{ -Download output file from the API -} diff --git a/man/format_annotate_args.Rd b/man/format_annotate_args.Rd deleted file mode 100644 index 3c2a09b..0000000 --- a/man/format_annotate_args.Rd +++ /dev/null @@ -1,23 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/argument_parser.R -\name{format_annotate_args} -\alias{format_annotate_args} -\title{Format parsed arguments for ANNotate function} -\usage{ -format_annotate_args(args) -} -\arguments{ -\item{args}{List of parsed arguments from parse_annotate_args()} -} -\value{ -Named list of formatted arguments for ANNotate function -} -\description{ -Format parsed arguments for ANNotate function -} -\examples{ -\dontrun{ -args <- parse_annotate_args() -formatted_args <- format_annotate_args(args) -} -} diff --git a/man/get_data.Rd b/man/get_data.Rd deleted file mode 100644 index 8c3a629..0000000 --- a/man/get_data.Rd +++ /dev/null @@ -1,21 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/seurat_utils.R -\name{get_data} -\alias{get_data} -\title{Get data from a Seurat object layer} -\usage{ -get_data(object, assay, layer = "data") -} -\arguments{ -\item{object}{Seurat object} - -\item{assay}{Assay name} - -\item{layer}{Layer name (default: "data")} -} -\value{ -Layer data -} -\description{ -Get data from a Seurat object layer -} diff --git a/man/if_gpu.Rd b/man/if_gpu.Rd deleted file mode 100644 index 872e763..0000000 --- a/man/if_gpu.Rd +++ /dev/null @@ -1,14 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/gpu_utils.R -\name{if_gpu} -\alias{if_gpu} -\title{Check for NVIDIA GPU availability} -\usage{ -if_gpu() -} -\value{ -Logical indicating whether an NVIDIA GPU is available -} -\description{ -Check for NVIDIA GPU availability -} diff --git a/man/listen_to_progress.Rd b/man/listen_to_progress.Rd deleted file mode 100644 index 17b1d4e..0000000 --- a/man/listen_to_progress.Rd +++ /dev/null @@ -1,21 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/api_interface.R -\name{listen_to_progress} -\alias{listen_to_progress} -\title{Listen to progress updates from the API} -\usage{ -listen_to_progress(url, file_path, ...) -} -\arguments{ -\item{url}{API endpoint URL} - -\item{file_path}{Path to the file being processed} - -\item{...}{Additional arguments passed to the API} -} -\value{ -Logical indicating success (TRUE) or failure (FALSE) -} -\description{ -Listen to progress updates from the API -} diff --git a/man/package_obj.Rd b/man/package_obj.Rd deleted file mode 100644 index 5b9a77e..0000000 --- a/man/package_obj.Rd +++ /dev/null @@ -1,34 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/seurat_utils.R -\name{package_obj} -\alias{package_obj} -\title{Package embeddings into a Seurat object} -\usage{ -package_obj( - extract_embeddings, - embeddings_dict, - umap_embeddings, - umap_embeddings_dict, - query_cells_df, - query_obj -) -} -\arguments{ -\item{extract_embeddings}{Whether embeddings are extracted} - -\item{embeddings_dict}{Dictionary of embeddings} - -\item{umap_embeddings}{Whether UMAP embeddings are computed} - -\item{umap_embeddings_dict}{Dictionary of UMAP embeddings} - -\item{query_cells_df}{Cell metadata} - -\item{query_obj}{Seurat object} -} -\value{ -Updated Seurat object -} -\description{ -Package embeddings into a Seurat object -} diff --git a/man/parse_annotate_args.Rd b/man/parse_annotate_args.Rd deleted file mode 100644 index 7dddd59..0000000 --- a/man/parse_annotate_args.Rd +++ /dev/null @@ -1,19 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/argument_parser.R -\name{parse_annotate_args} -\alias{parse_annotate_args} -\title{Parse command line arguments for the ANNotate function} -\usage{ -parse_annotate_args() -} -\value{ -A list of parsed arguments -} -\description{ -Parse command line arguments for the ANNotate function -} -\examples{ -\dontrun{ -args <- parse_annotate_args() -} -} diff --git a/man/process_azimuth_results.Rd b/man/process_azimuth_results.Rd deleted file mode 100644 index 94ea356..0000000 --- a/man/process_azimuth_results.Rd +++ /dev/null @@ -1,17 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/api_interface.R -\name{process_azimuth_results} -\alias{process_azimuth_results} -\title{Process results from PanAzimuth API} -\usage{ -process_azimuth_results(results) -} -\arguments{ -\item{results}{Raw results from the API} -} -\value{ -Processed results in Seurat object format -} -\description{ -Process results from PanAzimuth API -} diff --git a/man/process_rds_file.Rd b/man/process_rds_file.Rd deleted file mode 100644 index d2cc9bc..0000000 --- a/man/process_rds_file.Rd +++ /dev/null @@ -1,18 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/cloud.R -\name{process_rds_file} -\alias{process_rds_file} -\title{Process RDS file through the cloud API} -\usage{ -process_rds_file(api_base_url, file_path, ...) -} -\arguments{ -\item{api_base_url}{Base URL for the API} - -\item{file_path}{Path to input RDS file} - -\item{...}{Additional arguments passed to the API} -} -\description{ -Process RDS file through the cloud API -} diff --git a/man/read_obj_min.Rd b/man/read_obj_min.Rd deleted file mode 100644 index c029a27..0000000 --- a/man/read_obj_min.Rd +++ /dev/null @@ -1,21 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/seurat_utils.R -\name{read_obj_min} -\alias{read_obj_min} -\title{Read and process a Seurat object} -\usage{ -read_obj_min(query_obj, feature_names_col, assay_default = "RNA") -} -\arguments{ -\item{query_obj}{Seurat object} - -\item{feature_names_col}{Column name for feature names} - -\item{assay_default}{Default assay name} -} -\value{ -List containing processed data -} -\description{ -Read and process a Seurat object -} diff --git a/man/read_seurat_object.Rd b/man/read_seurat_object.Rd deleted file mode 100644 index 68d8a85..0000000 --- a/man/read_seurat_object.Rd +++ /dev/null @@ -1,25 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/read_seurat_object.R -\name{read_seurat_object} -\alias{read_seurat_object} -\title{Read a Seurat object from an RDS file} -\usage{ -read_seurat_object(filepath, assay_name = "RNA") -} -\arguments{ -\item{filepath}{Character string specifying the path to the RDS file} - -\item{assay_name}{Character string specifying the assay to use (default: "RNA")} -} -\value{ -A valid Seurat object -} -\description{ -Reads a Seurat object from an RDS file, ensuring -compatibility with Seurat versions >= 4.4.0 -} -\examples{ -\dontrun{ -seu_obj <- read_seurat_object("path/to/seurat_object.rds") -} -} diff --git a/man/run_azimuth_api.Rd b/man/run_azimuth_api.Rd deleted file mode 100644 index 0fe8875..0000000 --- a/man/run_azimuth_api.Rd +++ /dev/null @@ -1,19 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/api_interface.R -\name{run_azimuth_api} -\alias{run_azimuth_api} -\title{Run the PanAzimuth API analysis} -\usage{ -run_azimuth_api(input_data, api_endpoint) -} -\arguments{ -\item{input_data}{Input data for analysis} - -\item{api_endpoint}{API endpoint URL} -} -\value{ -Analysis results -} -\description{ -Run the PanAzimuth API analysis -} diff --git a/man/safe_progress_stream.Rd b/man/safe_progress_stream.Rd deleted file mode 100644 index 1df5a8e..0000000 --- a/man/safe_progress_stream.Rd +++ /dev/null @@ -1,21 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/api_interface.R -\name{safe_progress_stream} -\alias{safe_progress_stream} -\title{Helper function for \code{listen_to_progress} -Safely process SSE stream with error handling} -\usage{ -safe_progress_stream(url, handle) -} -\arguments{ -\item{url}{API endpoint URL} - -\item{handle}{Curl handle with the appropriate form data set} -} -\value{ -Logical indicating success (TRUE) or failure (FALSE) -} -\description{ -Helper function for \code{listen_to_progress} -Safely process SSE stream with error handling -} diff --git a/man/save_seurat_object.Rd b/man/save_seurat_object.Rd deleted file mode 100644 index 640aa1e..0000000 --- a/man/save_seurat_object.Rd +++ /dev/null @@ -1,19 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/save_seurat_object.R -\name{save_seurat_object} -\alias{save_seurat_object} -\title{Save a Seurat object to an RDS file} -\usage{ -save_seurat_object(seurat_obj, filepath) -} -\arguments{ -\item{seurat_obj}{A Seurat object to save} - -\item{filepath}{Character string specifying the output path} -} -\value{ -Invisible filepath where the object was saved -} -\description{ -Saves a Seurat object to an RDS file after applying "keto diet" to reduce size -} diff --git a/tests/testthat.R b/tests/testthat.R new file mode 100644 index 0000000..7ef857c --- /dev/null +++ b/tests/testthat.R @@ -0,0 +1,4 @@ +library(testthat) +library(AzimuthAPI) + +test_check("AzimuthAPI") diff --git a/tests/testthat/test_interface.R b/tests/testthat/test_interface.R new file mode 100644 index 0000000..771942f --- /dev/null +++ b/tests/testthat/test_interface.R @@ -0,0 +1,66 @@ +# Note - these tests are designed to test interface and packaging, i.e. they do not interact with the real server to test annotation results + +make_test_object <- function() { + obj <- readRDS(test_path("test_obj.rds")) + return(Seurat::NormalizeData(obj, verbose = FALSE)) +} + +test_that("CloudAzimuth returns an annotated Seurat object", { + query <- make_test_object() + expected <- query + expected@meta.data$final_level_labels <- rep(c("T cell", "B cell"), length.out = ncol(expected)) + expected@meta.data$azimuth_label <- expected@meta.data$final_level_labels + + # dummy functions to avoid server calls + testthat::local_mocked_bindings( + check_api_version = function(api_base_url) { + invisible(api_base_url) + }, + process_rds_file = function(api_base_url, file_path, ...) { # + saveRDS(expected, file = sub("\\.rds$", "_ANN.rds", file_path)) + invisible(NULL) + }, + .package = "AzimuthAPI" + ) + + annotated <- CloudAzimuth(query, ip = "azimuthapi.satijalab.org", port = 5000) + + expect_s4_class(annotated, "Seurat") + expect_true(all(c("final_level_labels", "azimuth_label") %in% colnames(annotated@meta.data))) + expect_identical(unname(annotated@meta.data$azimuth_label), rep(c("T cell", "B cell"), length.out = ncol(query))) +}) + +test_that("ANNotate returns an annotated Seurat object", { + skip_if_not_installed("Seurat") + skip_if_not_installed("reticulate") + + query <- make_test_object() + + # dummy functions to avoid reticulate setup + fake_import <- function(module, ...) { + if (identical(module, "panhumanpy.ANNotate")) { + return(list( + annotate_core = function(X_query, query_features, cells_meta, annotation_pipeline, eval_batch_size, normalization_override, norm_check_batch_size, + output_mode, refine_labels, map_to_cl, include_cl_id, extract_embeddings, umap_embeddings, n_neighbors, n_components, + metric, min_dist, umap_lr, umap_seed, spread, verbose, init, model_version) { + cells_meta$final_level_labels <- rep(c("T cell", "B cell"), length.out = nrow(cells_meta)) + list(embeddings_dict = list(), umap_dict = list(), cells_meta = cells_meta) + } + )) + } + + if (identical(module, "scipy.sparse")) { + return(list(csr_matrix = function(x) x)) + } + + stop(sprintf("Unexpected module request: %s", module)) + } + + testthat::local_mocked_bindings(import = fake_import, r_to_py = function(x) x, .package = "reticulate") + + annotated <- ANNotate(query, extract_embeddings = FALSE, umap_embeddings = FALSE) + + expect_s4_class(annotated, "Seurat") + expect_true(all(c("final_level_labels", "azimuth_label") %in% colnames(annotated@meta.data))) + expect_identical(unname(annotated@meta.data$azimuth_label), rep(c("T cell", "B cell"), length.out = ncol(query))) +}) diff --git a/test_obj.rds b/tests/testthat/test_obj.rds similarity index 100% rename from test_obj.rds rename to tests/testthat/test_obj.rds diff --git a/vignettes/intro_azimuthapi_vignette.Rmd b/vignettes/intro_azimuthapi_vignette.Rmd new file mode 100644 index 0000000..d3be531 --- /dev/null +++ b/vignettes/intro_azimuthapi_vignette.Rmd @@ -0,0 +1,238 @@ +--- +title: "Pan-human Azimuth in R" +output: + html_document: + toc: true + toc_depth: 3 + theme: united + df_print: kable +date: 'Compiled: `r format(Sys.Date(), "%B %d, %Y")`' +--- + +```{r setup, include=FALSE} +all_times <- list() # store the time for each chunk +knitr::knit_hooks$set(time_it = local({ + now <- NULL + function(before, options) { + if (before) { + now <<- Sys.time() + } else { + res <- difftime(Sys.time(), now, units = "secs") + all_times[[options$label]] <<- res + } + } +})) + +knitr::opts_chunk$set( + tidy = TRUE, + tidy.opts = list(width.cutoff = 95), + fig.align = "center", + message = FALSE, + warning = FALSE, + time_it = TRUE, + error = TRUE +) +``` + +This vignette demonstrates how to install and use `AzimuthAPI`, an interface for cell type annotation using the Pan-human Azimuth neural network. + +The API takes a Seurat object and returns hierarchical cell type predictions, confidence scores, and low-dimensional embeddings useful for visualization. + +# Install AzimuthAPI + +The Pan-human Azimuth R API, `AzimuthAPI`, is available for installation via Github: + +```{r install.panazimuth, eval=FALSE} +remotes::install_github("satijalab/AzimuthAPI") +``` + +# Load packages & data + +```{r load.packages} +library(AzimuthAPI) +library(Seurat) +library(SeuratData) +library(ggplot2) +library(dplyr) +``` + +For this vignette, we will demonstrate annotating a dataset of human bone marrow mononuclear (BMNC) cells that we published as part of [Stuart, Butler et al Cell 2019](https://doi.org/10.1016/j.cell.2019.05.031). This dataset is available through `SeuratData`. + +```{r install.data, eval = FALSE} +bmcite <- InstallData("bmcite") +``` + +Load the data from `SeuratData`. Here we also make sure to normalize the data before running Pan-human Azimuth. + +```{r load.data} +bmcite <- LoadData("bmcite") +bmcite <- NormalizeData(bmcite) +``` + +# Annotate a dataset with Pan-human Azimuth + +Users have two options for running Pan-human Azimuth on their dataset: cloud-based or local. The cloud-based option is the easiest to use, and requires no additional setup. The local option requires setting up a Python environment with panhumanpy (Python package for the Pan-human Azimuth neural network) and its dependencies. + +Annotation results, regardless of the method used, are stored in the Seurat object as cell-level metadata, and the embeddings generated by the underlying neural network are stored in a new Seurat reduction. + +## Option 1: run `CloudAzimuth` + +The `CloudAzimuth` function runs Pan-human Azimuth cell type annotation on a Seurat object via a cloud-based API deployed on AWS. + +```{r run.panazimuth, results = 'hide', message = FALSE, warning = FALSE} +bmcite_ann <- CloudAzimuth(bmcite) +``` + +## Option 2: run `ANNotate` + +First, set up a local Python environment: + +```{bash setup.local.env, eval = FALSE} +# Create a new conda environment for Pan-human Azimuth +conda create -n pan-human-azimuth python=3.9 +# Activate the environment +conda activate pan-human-azimuth +# Install Pan-human Azimuth and its dependencies +# for more installation options, see https://satijalab.org/pan_human_azimuth/python/panhumanpy.html +pip install panhumanpy +``` + +Then ask `reticulate` to use that environment: + +```{r use.local.env} +library(reticulate) +# Use the local Python environment for Pan-human Azimuth +use_condaenv("pan-human-azimuth", required = TRUE) +``` + +Finally, run the `ANNotate` function to annotate your dataset. Note that the initial run may take longer, as the model weights are downloaded and cached locally. + +```{r run.local.panazimuth, results = 'hide', message = FALSE, warning = FALSE} +bmcite_ann <- ANNotate(bmcite) +``` + +# Examine model output + +Pan-human Azimuth returns predictions in multiple formats, as well as softmax probability scores to estimate model confidence. + +For more information on the output, see the reference. + +```{r} +# Inspect cell-level metadata returned by the model +md <- bmcite_ann[[]] +# New columns - ex. final_level_labels, final_level_confidence, full_hierarchical_labels - have been added to the Seurat object metadata +colnames(md) +``` + +We can view a histogram of softmax probabilities, which reflect the model confidence associated with each cell's annotation: + +```{r softmax.histogram} +ggplot(md, aes(x = final_level_confidence)) + + geom_histogram(bins = 20, fill = "skyblue") + + labs(x = "Softmax probability", y = "Count", title = "Histogram of Softmax Probabilities") + + theme_bw() +``` + +As a QC / filtering step, we can remove cells with low scores: + +```{r subset.low.scores} +# remove cells with low scores +bmcite_qc <- subset(bmcite_ann, final_level_confidence > 0.5) +``` + +# Generate UMAP of Azimuth embeddings + +Pan-human Azimuth provides embeddings of the encoding layer generated by the underlying neural network, stored in the `azimuth_embed` 128-dimensional reduction. We can use this to generate a two-dimensional visualization of the dataset, to examine cell type predictions in this low-dimensional space. + +```{r run.umap} +# Use azimuth_embed as input to generate UMAP +bmcite_qc <- RunUMAP(bmcite_qc, dims = 1:128, reduction = "azimuth_embed", reduction.name = "azimuth_umap") +``` + +# Visualize predictions + +The full_hierarchical_label for each cell provides the model's classification at each level of granularity, with different levels separated by the '|' character. + +```{r pha.labels.1, fig.width = 6, fig.height = 6} +p2 <- DimPlot(bmcite_qc, group.by = "full_hierarchical_labels", label.size = 1.5, label = T, reduction = "azimuth_umap", + repel = TRUE) + NoLegend() +p2 +``` + +To avoid the long length of the full hierarchical label, we also output the final level of granularity in a separate column, `final_level_labels`, which is the last (most specific)level of the hierarchy for each cell. + +```{r pha.labels.2, fig.width = 6, fig.height = 6} +p3 <- DimPlot(bmcite_qc, group.by = "final_level_labels", label.size = 2, label = T, reduction = "azimuth_umap") + + NoLegend() +p3 +``` + +```{r pha.confidence, fig.width = 6, fig.height = 6} +p4 <- FeaturePlot(bmcite_qc, features = "final_level_confidence", reduction = "azimuth_umap") +p4 +``` + +We also postprocess our predictions to provide labels at three consistent levels of granularity for easy handling, marking any cell with an invalid full hierarchical label (based on `full_consistent_hierarchy`) as `False`. + +- `azimuth_broad`: Corresponds to level_zero_labels (e.g. Immune cell) +- `azimuth_medium`: Medium level of granularity (e.g. T cell) +- `azimuth_fine`: High level of granularity (e.g. Treg cell) + +These categories provide a consistent level of granularity for each cell, but may differ from the `final_level_label`, either by forcing the model to predict further along the cell type hierarchy than its intial prediction, or by rolling back its prediction to a lower level of granularity. + +```{r refined.labels.medium, fig.width = 6, fig.height = 6} +p5 <- DimPlot(bmcite_qc, group.by = "azimuth_medium", label.size = 3, label = T, reduction = "azimuth_umap") + + NoLegend() +p5 +``` + +```{r refined.labels.fine, fig.width = 6, fig.height = 6} +p6 <- DimPlot(bmcite_qc, group.by = "azimuth_fine", label.size = 3, label = T, reduction = "azimuth_umap") + + NoLegend() +p6 +``` + +To remove the number of labels displayed, we can filter labels with less than a certain number of cells using `PrepLabel`. + +Here, we filter labels with less than 20 cells per label. This can be useful to filter outliers, especially as Pan-human Azimuth does not perform smoothing of single-cell labels by cluster. Therefore, a single outlier annotation for one cell will still display the outlier label on a visualization, and the `PrepLabel` function can help with this. + +```{r prepped.labels, fig.width = 6, fig.height = 6} +bmcite_qc <- PrepLabel(bmcite_qc, "azimuth_fine", "azimuth_fine_filtered", cutoff = 20) +p7 <- DimPlot(bmcite_qc, group.by = "azimuth_fine_filtered", label.size = 3, label = T, reduction = "azimuth_umap") + + NoLegend() +p7 +``` + +# Visualize differentially expressed features + +The `make_azimuth_QC_heatmaps` function allows you to easily explore the quality of predicted labels by creating expression heatmaps by predicted cell type, with optional parameters for improved visualization. + +Plots are saved by `azimuth_broad` categories by default, with the exception of immune cell types grouped separately by lymphoid or myeloid/erythroid subpopulations. + +```{r make.heatmaps} +plots <- make_azimuth_QC_heatmaps(bmcite_qc) +print(length(plots)) +``` + +```{r lymphoid.heatmap, fig.height = 12, fig.width = 12} +p8 <- plots[['Immune_Lymphoid cell_1']] +print(p8) +``` + +```{r myeloid.heatmap, fig.height = 12, fig.width = 12} +p9 <- plots[['Immune_Myeloid cell_1']] +print(p9) +``` + +# More information + +- For detailed options and usage of the `AzimuthAPI` functions, see the AzimuthAPI [function reference](https://google.com). +- For another example of the usage of this package, see the [Analysis and visualization of Visium HD data with cell segmentations](https://satijalab.org/seurat/articles/visiumhd_analysis_cell_segmentations) vignette on the Seurat website. + +
**Session Info** + +```{r session.info} +sessionInfo() +``` + +
\ No newline at end of file diff --git a/vignettes/quickstart_vignette.Rmd b/vignettes/quickstart_vignette.Rmd deleted file mode 100644 index 7f36563..0000000 --- a/vignettes/quickstart_vignette.Rmd +++ /dev/null @@ -1,147 +0,0 @@ ---- -title: "Pan-human Azimuth R API Vignette" -output: - html_document: - theme: united - df_print: kable -date: 'Compiled: `r format(Sys.Date(), "%B %d, %Y")`' ---- -```{r setup, include=FALSE} -all_times <- list() # store the time for each chunk -knitr::knit_hooks$set(time_it = local({ - now <- NULL - function(before, options) { - if (before) { - now <<- Sys.time() - } else { - res <- difftime(Sys.time(), now, units = "secs") - all_times[[options$label]] <<- res - } - } -})) -knitr::opts_chunk$set( - tidy = TRUE, - tidy.opts = list(width.cutoff = 95), - fig.align = "center", - message = FALSE, - warning = FALSE, - time_it = TRUE, - error = TRUE -) -``` - -This vignette demonstrates how to install and run the Pan-human Azimuth R Cloud API. The API takes a Seurat object and returns hierarchical cell type predictions, confidence scores, and a low-dimensional visualization. - -# Setup the Environment -Install the Pan-human Azimuth R API from Github -```{r install.panazimuth, eval=FALSE} -devtools::install_github("satijalab/AzimuthAPI", auth_token = "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") -``` - -# Run ANNotate -For this vignette, we will demonstrate annotating a dataset of human bone marrow mononuclear (BMNC) cells that we published as part of (Stuart, Butler et al Cell 2019)[https://www.cell.com/cell/fulltext/S0092-8674(19)30559-8]. This dataset is available through SeuratData. -```{r load.packages} -library(Seurat) -library(SeuratData) -library(AzimuthAPI) -library(dplyr) -``` -```{r install.data, eval = FALSE} -bmcite <- InstallData("bmcite") -``` - -```{r load.data} -bmcite <- LoadData("bmcite") -``` - -The `CloudAzimuth` function runs Pan-human Azimuth cell type predictions on a Seurat object via a cloud-based API deployed on AWS, and returns results stored in the object's cell-level metadata. - -```{r run.panazimuth, results = FALSE, cache=TRUE} -bmcite <- CloudAzimuth(bmcite) -``` - -# Visualize predictions - -You can visualize the predicted cell type labels on `azimuthumap`, a UMAP visualization of the neural network model's encoding layer. - -Pan-human Azimuth returns predictions in multiple formats, as well as softmax probability scores to estimate model confidence: - -* `full_hierarchical_labels`: The full predicted label for each cell, where '|' delimits each hierarchical level -* `level_zero_labels`: The most broad level of the full hierarchical label -* `final_level_labels`: The most granular level of the full hierarchical label -* `final_level_softmax_p*rob`: The model’s predicted probability for each cell’s assigned `final_level_label`, ranging from 0 to 1 -* `full_consistent_hierarchy`: A flag that indicates whether all levels of the hierarchy are internally consistent - -```{r} -# All cell-level metadata that is returned -# Note that the first 10 columns existed in the object prior to running Azimuth -colnames(bmcite[[]]) -``` -```{r, include = FALSE} -#bmcite <- PrepLabel(bmcite, "full_hierarchical_labels", "full_hierarchical_labels", cutoff=20) -#bmcite <- PrepLabel(bmcite, "final_level_labels", "final_level_labels", cutoff=20) -``` - -The `full_hierarchical_label` for each cell provides the model's classification at each level of granularity, with different levels separated by the '|' character. - -```{r panazimuth.labels.1} -p2 <- DimPlot(bmcite, group.by = "full_hierarchical_labels", label.size = 1.5, repel=T, label=T) + NoLegend() -p2 -``` -To avoid the long length of the full hierarchical label, we also output the final level of granularity - -```{r panazimuth.labels.2, fig.width = 9} -p3 <- DimPlot(bmcite, group.by = "final_level_labels", label.size = 3, repel=T, label=T) + NoLegend() -p4 <- FeaturePlot(bmcite, features = "final_level_softmax_prob") -p3 | p4 -``` - -We also postprocess our predictions to provide labels at three consistent levels of granularity for easy handling, marking any cell with an invalid full hierarchical label (based on `full_consistent_hierarchy`) as `False`. - -* `azimuth_broad`: Corresponds to level_zero_labels (i.e. Immune cell) -* `azimuth_medium`: Medium level of granularity (i.e. T cell) -* `azimuth_fine`: High level of granularity (i.e. Treg cell) - -These categories provide a consistent level of granularity for each cell, but may differ from the `final_level_label`, either by forcing the model to predict further along the cell type hierarchy than its intial prediction, or by rolling back its prediction to a lower level of granularity. - -```{r refined.labels.medium} -p5 <- DimPlot(bmcite, group.by = "azimuth_medium", label.size = 3, repel=T, label=T) + NoLegend() -p5 -``` -```{r refined.labels.fine} -p6 <- DimPlot(bmcite, group.by = "azimuth_fine", label.size = 3, repel=T, label=T) + NoLegend() -p6 -``` - -To remove the number of labels displayed, you can filter labels with less than a certain number of cells using `PrepLabel`. - -```{r prepped.labels} -bmcite <- PrepLabel(bmcite, "azimuth_fine", "azimuth_fine_filtered", cutoff=20) -p7 <- DimPlot(bmcite, group.by = "azimuth_fine_filtered", label.size = 3, repel=T, label=T) + NoLegend() -p7 -``` - -# Visualize differentially expressed features -The `make_azimuth_QC_heatmaps` function allows you to easily explore the quality of predicted labels by creating expression heatmaps by predicted cell type, with optional parameters for improved visualization: - -* `final_name`: Name of metadata columns to group cells by (default is `azimuth_fine`) -* `min.final.group`: Minimum number of cells under a cell type to be displayed -* `max.ids.per.plot`: Number of cell type labels displayed per plot -* `reorder`: Flag to indicate whether to reorder cell types by transcriptional similarity -* `cells.order`: Cell names to specify order of cells (i.e. by softmax probability) -* `save_folder_path`: Save plots as PNG files under specified folder path - -Plots are saved by `azimuth_broad` categories by default, with the exception of immune cell types grouped separately by lymphoid or myeloid/erythroid subpopulations. - -```{r make.heatmaps} -plots <- make_azimuth_QC_heatmaps(bmcite) -print(length(plots)) -``` -```{r lymphoid.heatmap, fig.height = 12,fig.width=12} -p8 <- plots[['Immune_Lymphoid cell_1']] -print(p8) -``` -```{r myeloid.heatmap, fig.height = 12,fig.width=12} -p9 <- plots[['Immune_Myeloid cell_1']] -print(p9) -``` \ No newline at end of file