From df21024eed7a1d4d98b839145734bd02c67d13eb Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Wed, 7 May 2025 16:48:54 +0300 Subject: [PATCH 01/28] init --- R/AllGenerics.R | 5 ++++ R/plotOrdination.R | 61 ++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 66 insertions(+) create mode 100644 R/plotOrdination.R diff --git a/R/AllGenerics.R b/R/AllGenerics.R index b14d864a..87d2e071 100644 --- a/R/AllGenerics.R +++ b/R/AllGenerics.R @@ -150,3 +150,8 @@ setGeneric("plotBarplot", signature = c("x"), function(x, ...) #' @export setGeneric("plotBoxplot", signature = c("object"), function(object, ...) standardGeneric("plotBoxplot")) + +#' @rdname plotOrdination +#' @export +setGeneric("plotOrdination", signature = c("x"), function(x, ...) + standardGeneric("plotOrdination")) diff --git a/R/plotOrdination.R b/R/plotOrdination.R new file mode 100644 index 00000000..8e98a5c1 --- /dev/null +++ b/R/plotOrdination.R @@ -0,0 +1,61 @@ +#' @name +#' plotOrdination +#' +#' @title +#' Create ordination plot +#' +#' @description +#' Ordinaton plotter +#' +#' @details +#' Creates ordination plot +#' +#' @return +#' A \code{ggplot2} object. +#' +#' @param x a +#' \code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}} +#' object. +#' +#' @param ... Additional parameters for plotting. +#' \itemize{ +#' \item \code{colour.by}: \code{NULL} or \code{character scalar}. Specifies a +#' variable from \code{colData(x)} or \code{rowData(x)} which is used to +#' colour observations. (Default: \code{NULL}) +#' } +#' +#' @examples +#' data("Tito2024QMP") +#' tse <- Tito2024QMP +#' +#' @seealso +#' \itemize{ +#' \item \code{\link[scater:plotReducedDim]{scater::plotReducedDim}} +#' } +#' +NULL + +#' @rdname plotOrdination +#' @export +setMethod("plotOrdination", signature = c(x = "SummarizedExperiment"), + function(x, ...){ + temp <- .check_ordination_input(x, ...) + df <- .get_ordination_data(x, ...) + p <- .ordination_plotter(df, ...) + return(p) + } +) + +################################ HELP FUNCTIONS ################################ + +.check_ordination_input <- function(x, ...){ + +} + +.get_ordination_data <- function(x, ...){ + +} + +.ordination_plotter <- function(df, ...){ + +} From 442d73dc8ec8650639947724f3ff81ea0d1ee362 Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Thu, 8 May 2025 09:38:32 +0300 Subject: [PATCH 02/28] up --- R/plotOrdination.R | 54 +++++++++++++++++++++++++++++++++++++++++----- R/utils.R | 9 ++++++-- 2 files changed, 56 insertions(+), 7 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 8e98a5c1..66dfa022 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -37,8 +37,8 @@ NULL #' @rdname plotOrdination #' @export -setMethod("plotOrdination", signature = c(x = "SummarizedExperiment"), - function(x, ...){ +setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), + function(x, dimred, colour.by = NULL, ...){ temp <- .check_ordination_input(x, ...) df <- .get_ordination_data(x, ...) p <- .ordination_plotter(df, ...) @@ -48,12 +48,56 @@ setMethod("plotOrdination", signature = c(x = "SummarizedExperiment"), ################################ HELP FUNCTIONS ################################ -.check_ordination_input <- function(x, ...){ +.check_ordination_input <- function( + x, dimred, + ncomponents = 2L, + colour.by = color.by, color.by = NULL, + shape.by = NULL, + size.by = NULL, + group.by = NULL, + pair.by = NULL, + order.by = NULL, + assay.type = "counts", + ...){ + # Check if there are any reduced dim present + if( lenght(reducedDims(x)) == 0L ){ + stop("No data present in reducedDim(x).", call. = FALSE) + } + # Check that dimred can be found + is_name <- .is_a_string(dimred) && dimred %in% reducedDimNames(x) + is_index <- .is_an_integer(dimred) && dimred > 0L && dimred <= length(reducedDims(x)) + if( !(is_name || is_index) ){ + stop("'dimred' must specify data from reducedDim(x). It must be one ", + "of the following options: '", + paste0(reducedDimNames(x), collapse = "', '"), "'", call. = FALSE) + } + # Check that ncomponents is correct. We can only visualize 2 components. + if( .is_an_integer(ncomponents) ){ + ncomponents <- seq_len(ncomponents) + } + if( !(.is_integer(ncomponents) && length(ncomponents) == 2L) ){ + stop("'ncomponents' must specify intger values.", call. = FALSE) + } -} + # Check aesthetic variables + temp <- .check_metadata_variable(tse, colour.by, TRUE, TRUE, FALSE, TRUE) + temp <- .check_metadata_variable(tse, shape.by, TRUE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, size.by, TRUE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, group.by, TRUE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, pair.by, TRUE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, order.by, TRUE, TRUE, FALSE, FALSE) + # If colour.by specifies rowname, we check assay.type as the abundance + # values are used for coloring + if( colour.by &in& rownames(x) ){ + temp <- .check_assay_present(assay.type, x) + } + # Check other flags -.get_ordination_data <- function(x, ...){ +} +.get_ordination_data <- function(x, dimred, ...){ + df <- reducedDim(x, dimred) |> as.data.frame() + df <- df[, ncomponents] } .ordination_plotter <- function(df, ...){ diff --git a/R/utils.R b/R/utils.R index 33daeaec..75ee2456 100644 --- a/R/utils.R +++ b/R/utils.R @@ -17,6 +17,7 @@ .is_function <- mia:::.is_function .get_name_in_parent <- mia:::.get_name_in_parent .is_an_integer <- mia:::.is_an_integer +.is_integer <- mia:::.is_integer TAXONOMY_RANKS <- mia:::TAXONOMY_RANKS .is_a_numeric <- mia:::.is_a_numeric .capitalize <- mia:::.capitalize @@ -55,8 +56,9 @@ TAXONOMY_RANKS <- mia:::TAXONOMY_RANKS } # This function checks whether variable can be found from colData or rowData. +# Optionally, we can look also from rownames. .check_metadata_variable <- function( - tse, var, row = FALSE, col = FALSE, multiple = FALSE, + tse, var, row = FALSE, col = FALSE, multiple = FALSE, rownames = FALSE, var.name = .get_name_in_parent(var)){ if( !.is_a_bool(multiple) ){ stop("'multiple' must be TRUE or FALSE.", call. = FALSE) @@ -68,11 +70,14 @@ TAXONOMY_RANKS <- mia:::TAXONOMY_RANKS check_values <- c() check_values <- c(check_values, if(col) colnames(colData(tse))) check_values <- c(check_values, if(row) colnames(rowData(tse))) + msg_values <- check_values + # Optionally, we look also rownames + check_values <- c(check_values, if(rownames) rownames(tse)) var_found <- all( var %in% check_values ) if( !(is_string && var_found) ){ stop("'", var.name, "' must be ", ifelse(multiple, "", "a single "), "character value from the following options: '", - paste0(check_values, collapse = "', '"), "'", call. = FALSE) + paste0(msg_values, collapse = "', '"), "'", call. = FALSE) } } return(NULL) From 44f58039cb6bf3ac814b0ebf1b7c20d5d583704b Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Thu, 8 May 2025 10:57:28 +0300 Subject: [PATCH 03/28] up --- R/plotOrdination.R | 72 +++++++++++++++++++++++++++++++++++++--------- 1 file changed, 58 insertions(+), 14 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 66dfa022..cb0aed3a 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -39,8 +39,9 @@ NULL #' @export setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), function(x, dimred, colour.by = NULL, ...){ - temp <- .check_ordination_input(x, ...) - df <- .get_ordination_data(x, ...) + args <- .check_ordination_input(x, dimred, colour.by = colour.by, ...) + df <- do.call(.get_ordination_data, args) + return(df) p <- .ordination_plotter(df, ...) return(p) } @@ -60,7 +61,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), assay.type = "counts", ...){ # Check if there are any reduced dim present - if( lenght(reducedDims(x)) == 0L ){ + if( length(reducedDims(x)) == 0L ){ stop("No data present in reducedDim(x).", call. = FALSE) } # Check that dimred can be found @@ -75,29 +76,72 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( .is_an_integer(ncomponents) ){ ncomponents <- seq_len(ncomponents) } - if( !(.is_integer(ncomponents) && length(ncomponents) == 2L) ){ - stop("'ncomponents' must specify intger values.", call. = FALSE) + if( !(.is_integer(ncomponents) && length(ncomponents) == 2L && all(ncomponents > 0L & ncomponents <= ncol(reducedDim(x, dimred)))) ){ + stop("'ncomponents' must specify columns from reducedDim(x, dimred) with integer values.", call. = FALSE) } # Check aesthetic variables - temp <- .check_metadata_variable(tse, colour.by, TRUE, TRUE, FALSE, TRUE) - temp <- .check_metadata_variable(tse, shape.by, TRUE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, size.by, TRUE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, group.by, TRUE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, pair.by, TRUE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, order.by, TRUE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, colour.by, FALSE, TRUE, FALSE, TRUE) + temp <- .check_metadata_variable(tse, shape.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, size.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, group.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, pair.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, order.by, FALSE, TRUE, FALSE, FALSE) # If colour.by specifies rowname, we check assay.type as the abundance # values are used for coloring - if( colour.by &in& rownames(x) ){ + if( !is.null(colour.by) && colour.by %in% rownames(x) ){ temp <- .check_assay_present(assay.type, x) } # Check other flags + # Put all the arguments into list that be fed to the data retrieval function + args <- list( + x = x, dimred = dimred, + ncomponents = ncomponents, + colour.by = colour.by, + shape.by = shape.by, + size.by = size.by, + group.by = group.by, + pair.by = pair.by, + order.by = order.by, + assay.type = "counts" + ) + args <- c(args, list(...)) + return(args) } -.get_ordination_data <- function(x, dimred, ...){ - df <- reducedDim(x, dimred) |> as.data.frame() +.get_ordination_data <- function(x, dimred, ncomponents, colour.by, shape.by, size.by, group.by, pair.by, order.by, assay.type, ...){ + df <- reducedDim(x, dimred) + if( is.null(colnames(df)) ){ + colnames(df) <- paste0(dimred, seq_len(ncol(df))) + } + df <- df |> as.data.frame() df <- df[, ncomponents] + x_var <- colnames(df)[[1L]] + y_var <- colnames(df)[[2L]] + + cols <- c(shape.by, size.by, group.by, pair.by, order.by) + if( !is.null(colour.by) && colour.by %in% colnames(colData(x)) ){ + cols <- c(cols, colour.by) + } else if( !is.null(colour.by) ){ + df[[colour.by]] <- assay(x, assay.type)[colour.by, ] + } + cd <- colData(x)[, cols, drop = FALSE] + df <- cbind(df, cd) + + attributes(df) <- c( + attributes(df), + x = x_var, + y = y_var, + colour.by = colour.by, + shape.by = shape.by, + size.by = size.by, + group.by = group.by, + pair.by = pair.by, + order.by = order.by, + assay.type = assay.type + ) + return(df) } .ordination_plotter <- function(df, ...){ From 11916b463d6a2b6fe756f2a5a2aa1adebc445a6f Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Fri, 9 May 2025 09:39:50 +0300 Subject: [PATCH 04/28] up --- R/plotBoxplot.R | 7 +++--- R/plotOrdination.R | 56 +++++++++++++++++++++++++++++++++++++++++----- 2 files changed, 54 insertions(+), 9 deletions(-) diff --git a/R/plotBoxplot.R b/R/plotBoxplot.R index 2a38724a..74b4098a 100644 --- a/R/plotBoxplot.R +++ b/R/plotBoxplot.R @@ -710,14 +710,15 @@ setMethod("plotBoxplot", signature = c(object = "SummarizedExperiment"), # This function adds points to plot .add_points_layer <- function( - p, df, point.alpha = 0.65, point.size = 2, point.shape = 21, + p, df, x = "x_point", y = "y_point", + point.alpha = 0.65, point.size = 2, point.shape = 21, point.colour = point.color, point.color = "grey70", ...){ # To disable "no visible binding for global variable" message in cmdcheck x_point <- y_point <- NULL args <- list( mapping = aes( - x = x_point, - y = y_point, + x = .data[[x]], + y = .data[[y]], colour = if(!is.null(attributes(df)[["colour.by"]])) .data[[attributes(df)[["colour.by"]]]], shape = if(!is.null(attributes(df)[["shape.by"]])) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index cb0aed3a..8dd6f316 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -38,10 +38,9 @@ NULL #' @rdname plotOrdination #' @export setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), - function(x, dimred, colour.by = NULL, ...){ + function(x, dimred, colour.by = color.by, color.by = NULL, ...){ args <- .check_ordination_input(x, dimred, colour.by = colour.by, ...) df <- do.call(.get_ordination_data, args) - return(df) p <- .ordination_plotter(df, ...) return(p) } @@ -53,11 +52,13 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), x, dimred, ncomponents = 2L, colour.by = color.by, color.by = NULL, + fill.by = NULL, shape.by = NULL, size.by = NULL, group.by = NULL, pair.by = NULL, order.by = NULL, + facet.by = NULL, assay.type = "counts", ...){ # Check if there are any reduced dim present @@ -82,11 +83,13 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # Check aesthetic variables temp <- .check_metadata_variable(tse, colour.by, FALSE, TRUE, FALSE, TRUE) + temp <- .check_metadata_variable(tse, fill.by, FALSE, TRUE, FALSE, TRUE) temp <- .check_metadata_variable(tse, shape.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, size.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, group.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, pair.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, order.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, facet.by, FALSE, TRUE, FALSE, FALSE) # If colour.by specifies rowname, we check assay.type as the abundance # values are used for coloring if( !is.null(colour.by) && colour.by %in% rownames(x) ){ @@ -99,33 +102,41 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), x = x, dimred = dimred, ncomponents = ncomponents, colour.by = colour.by, + fill.by = fill.by, shape.by = shape.by, size.by = size.by, group.by = group.by, pair.by = pair.by, order.by = order.by, + facet.by = facet.by, assay.type = "counts" ) args <- c(args, list(...)) return(args) } -.get_ordination_data <- function(x, dimred, ncomponents, colour.by, shape.by, size.by, group.by, pair.by, order.by, assay.type, ...){ +.get_ordination_data <- function(x, dimred, ncomponents, colour.by, fill.by, shape.by, size.by, group.by, pair.by, order.by, facet.by, assay.type, ...){ df <- reducedDim(x, dimred) if( is.null(colnames(df)) ){ - colnames(df) <- paste0(dimred, seq_len(ncol(df))) + colnames(df) <- paste0(dimred, ncomponents) } df <- df |> as.data.frame() df <- df[, ncomponents] x_var <- colnames(df)[[1L]] y_var <- colnames(df)[[2L]] - cols <- c(shape.by, size.by, group.by, pair.by, order.by) + cols <- c(shape.by, size.by, group.by, pair.by, order.by, facet.by) if( !is.null(colour.by) && colour.by %in% colnames(colData(x)) ){ cols <- c(cols, colour.by) } else if( !is.null(colour.by) ){ df[[colour.by]] <- assay(x, assay.type)[colour.by, ] } + if( !is.null(fill.by) && fill.by %in% colnames(colData(x)) ){ + cols <- c(cols, fill.by) + } else if( !is.null(fill.by) ){ + df[[fill.by]] <- assay(x, assay.type)[fill.by, ] + } + cols <- cols |> unique() cd <- colData(x)[, cols, drop = FALSE] df <- cbind(df, cd) @@ -134,16 +145,49 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), x = x_var, y = y_var, colour.by = colour.by, + fill.by = fill.by, shape.by = shape.by, size.by = size.by, group.by = group.by, pair.by = pair.by, order.by = order.by, + facet.by = facet.by, assay.type = assay.type ) return(df) } -.ordination_plotter <- function(df, ...){ +.ordination_plotter <- function(df, scales = "fixed", ...){ + p <- ggplot(df, aes( + x = .data[[attributes(df)[["x"]]]], + y = .data[[attributes(df)[["y"]]]], + colour = if( !is.null(attributes(df)[["colour.by"]]) ) + .data[[attributes(df)[["colour.by"]]]] + )) + p <- .add_points_layer(p, df, x = attributes(df)[["x"]], y = attributes(df)[["y"]], ...) + + # If facetting was specified, split plot to separate panels + if( !is.null(attributes(df)[["facet.by"]]) ){ + p <- p + + facet_wrap( + ~ .data[[attributes(df)[["facet.by"]]]], + scales = scales + ) + } + p <- p + theme_classic() + p <- p + labs(x = attributes(df)[["x"]]) + if( !is.null(attributes(df)[["fill.by"]]) ){ + p <- p + labs(fill = attributes(df)[["fill.by"]]) + } + if( !is.null(attributes(df)[["colour.by"]]) ){ + p <- p + labs(colour = attributes(df)[["colour.by"]]) + } + if( !is.null(attributes(df)[["shape.by"]]) ){ + p <- p + labs(shape = attributes(df)[["shape.by"]]) + } + if( !is.null(attributes(df)[["size.by"]]) ){ + p <- p + labs(shape = attributes(df)[["size.by"]]) + } + return(p) } From 1871b1e1bc701ce1ac6e541a712d29e306995fad Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Mon, 12 May 2025 12:15:31 +0300 Subject: [PATCH 05/28] up --- R/plotOrdination.R | 355 +++++++++++++++++++++++++++++++++++++++++++-- 1 file changed, 342 insertions(+), 13 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 8dd6f316..ce3b3e34 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -48,6 +48,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), ################################ HELP FUNCTIONS ################################ +# This method checks that the input is in correct format .check_ordination_input <- function( x, dimred, ncomponents = 2L, @@ -57,9 +58,12 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), size.by = NULL, group.by = NULL, pair.by = NULL, - order.by = NULL, + sort.by = NULL, facet.by = NULL, assay.type = "counts", + add.points = TRUE, add.ellipse = FALSE, add.density = FALSE, + add.centroids = FALSE, add.centroids.lines = FALSE, add.vectors = FALSE, + add.rotation = FALSE, ...){ # Check if there are any reduced dim present if( length(reducedDims(x)) == 0L ){ @@ -67,7 +71,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Check that dimred can be found is_name <- .is_a_string(dimred) && dimred %in% reducedDimNames(x) - is_index <- .is_an_integer(dimred) && dimred > 0L && dimred <= length(reducedDims(x)) + is_index <- .is_an_integer(dimred) && dimred > 0L && + dimred <= length(reducedDims(x)) if( !(is_name || is_index) ){ stop("'dimred' must specify data from reducedDim(x). It must be one ", "of the following options: '", @@ -77,8 +82,11 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( .is_an_integer(ncomponents) ){ ncomponents <- seq_len(ncomponents) } - if( !(.is_integer(ncomponents) && length(ncomponents) == 2L && all(ncomponents > 0L & ncomponents <= ncol(reducedDim(x, dimred)))) ){ - stop("'ncomponents' must specify columns from reducedDim(x, dimred) with integer values.", call. = FALSE) + if( !(.is_integer(ncomponents) && length(ncomponents) == 2L && + all(ncomponents > 0L & ncomponents <= ncol(reducedDim(x, dimred)))) + ){ + stop("'ncomponents' must specify columns from reducedDim(x, dimred) ", + "with integer values.", call. = FALSE) } # Check aesthetic variables @@ -88,7 +96,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), temp <- .check_metadata_variable(tse, size.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, group.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, pair.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, order.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(tse, sort.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, facet.by, FALSE, TRUE, FALSE, FALSE) # If colour.by specifies rowname, we check assay.type as the abundance # values are used for coloring @@ -96,7 +104,34 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), temp <- .check_assay_present(assay.type, x) } # Check other flags + if( !.is_a_bool(add.points) ){ + stop("'add.points' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.ellipse) ){ + stop("'add.ellipse' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.density) ){ + stop("'add.density' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.centroids) ){ + stop("'add.centroids' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.centroids.lines) ){ + stop("'add.centroids.lines' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.vectors) ){ + stop("'add.vectors' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(add.rotation) ){ + stop("'add.rotation' must be TRUE or FALSE.", call. = FALSE) + } + # add.density cannot be specified simultaneously with fill.by as they both + # are using fill aesthetic and we can have only on fill scale. + if( add.density && !is.null(fill.by) ){ + stop("Both 'add.density' and 'fill.by' cannot be specified ", + "simultaneously.", call. = FALSE) + } # Put all the arguments into list that be fed to the data retrieval function args <- list( x = x, dimred = dimred, @@ -107,7 +142,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), size.by = size.by, group.by = group.by, pair.by = pair.by, - order.by = order.by, + sort.by = sort.by, facet.by = facet.by, assay.type = "counts" ) @@ -115,17 +150,42 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), return(args) } -.get_ordination_data <- function(x, dimred, ncomponents, colour.by, fill.by, shape.by, size.by, group.by, pair.by, order.by, facet.by, assay.type, ...){ +# This function retrieves the data from reducedDim +.get_ordination_data <- function( + x, dimred, ncomponents, colour.by, fill.by, shape.by, size.by, group.by, + pair.by, sort.by, facet.by, assay.type, ...){ + # Get data and store the original attributes that might include rotation + # data, for instance df <- reducedDim(x, dimred) + orig_attributes <- attributes(df) + orig_attributes <- orig_attributes[ + !names(orig_attributes) %in% c("dim", "dimnames") ] + + # Add colnames if they are not present if( is.null(colnames(df)) ){ colnames(df) <- paste0(dimred, ncomponents) } df <- df |> as.data.frame() + # Take only 2 specified columns df <- df[, ncomponents] x_var <- colnames(df)[[1L]] y_var <- colnames(df)[[2L]] - cols <- c(shape.by, size.by, group.by, pair.by, order.by, facet.by) + # Get rotation data adnd put it in correct format + rotation <- NULL + rotation_names <- c("rotation") + if( any(rotation_names %in% names(orig_attributes)) ){ + rotation_names <- rotation_names[[1L]] + rotation <- orig_attributes[[rotation_names]] |> as.data.frame() + rotation <- rotation[, ncomponents, drop = FALSE] + colnames(rotation) <- colnames(df) + } + + # List data that is fetched from colData + cols <- c(shape.by, size.by, group.by, pair.by, sort.by, facet.by) + # colour.by and fill.by can also specify a feature and its abundance. + # Either get their column name to fetch from colData or directly add them + # to data. if( !is.null(colour.by) && colour.by %in% colnames(colData(x)) ){ cols <- c(cols, colour.by) } else if( !is.null(colour.by) ){ @@ -136,10 +196,32 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } else if( !is.null(fill.by) ){ df[[fill.by]] <- assay(x, assay.type)[fill.by, ] } + # Get colData and merge it with the data cols <- cols |> unique() cd <- colData(x)[, cols, drop = FALSE] df <- cbind(df, cd) + # Sort the data. For instance, if we want to add arrows between consecutive + # time points, this is essential step. + if( !is.null(sort.by) ){ + df <- df[order(df[[sort.by]]), , drop = FALSE] + } + + # Calculate centroid data + grouping_var <- c( + colour.by, facet.by, group.by, fill.by) |> unique() + df_centroids <- df |> + group_by(across(all_of(grouping_var))) |> + summarise( + x_centroid = mean(.data[[x_var]], na.rm = TRUE), + y_centroid = mean(.data[[y_var]], na.rm = TRUE), + .groups = "drop" + ) + # Add global mean + df_centroids[["x_global"]] <- mean(df[[x_var]], na.rm = TRUE) + df_centroids[["y_global"]] <- mean(df[[y_var]], na.rm = TRUE) + + # Add additional information to attributes of df. attributes(df) <- c( attributes(df), x = x_var, @@ -150,22 +232,63 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), size.by = size.by, group.by = group.by, pair.by = pair.by, - order.by = order.by, + sort.by = sort.by, facet.by = facet.by, assay.type = assay.type ) + attr(df, "rotation") <- rotation + attr(df, "centroids") <- df_centroids return(df) } -.ordination_plotter <- function(df, scales = "fixed", ...){ +# This method is the main plotter function. +.ordination_plotter <- function( + df, scales = "fixed", add.points = TRUE, add.ellipse = FALSE, + add.density = FALSE, add.centroids = FALSE, add.centroids.lines = FALSE, + add.vectors = FALSE, add.rotation = FALSE, ...){ + # Initialize the plot p <- ggplot(df, aes( x = .data[[attributes(df)[["x"]]]], y = .data[[attributes(df)[["y"]]]], colour = if( !is.null(attributes(df)[["colour.by"]]) ) .data[[attributes(df)[["colour.by"]]]] )) - p <- .add_points_layer(p, df, x = attributes(df)[["x"]], y = attributes(df)[["y"]], ...) - + grouping_var <- c( + attributes(df)[["colour.by"]], attributes(df)[["facet.by"]], + attributes(df)[["group.by"]], attributes(df)[["fill.by"]]) |> unique() + # Add points, i.e., samples + if( add.points ){ + p <- .add_ordination_points(p, df, ...) + } + # Add ellipses + if( add.ellipse && length(grouping_var) > 0L ){ + p <- .add_ellipse_to_ordination(p, df, ...) + } + # Connect samples of subjects, groups etc + if( !is.null(attributes(df)[["pair.by"]]) ){ + p <- .connect_ordination_points(p, df) + } + # Add density to background + if( add.density ){ + # Background density + p <- .add_point_density(p, df) + } + # Add group centroids + if( add.centroids && length(grouping_var) > 0L ){ + p <- .add_centroids(p, df, ...) + } + # Add lines connecting points and group centroids + if( add.centroids.lines && length(grouping_var) > 0L ){ + p <- .add_centroid_lines(p, df, grouping_var, ...) + } + # Add species scores + if( add.rotation ){ + p <- .add_rotation(p, df) + } + # Add vectors from global centroid to group centroids + if( add.vectors && length(grouping_var) > 0L ){ + p <- .add_centroids_vector(p, df, ...) + } # If facetting was specified, split plot to separate panels if( !is.null(attributes(df)[["facet.by"]]) ){ p <- p + @@ -174,7 +297,213 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), scales = scales ) } + # Adjust theme + p <- .adjust_ordination_theme(p, df, ...) + return(p) +} +# Add points for ordination plot +.add_ordination_points <- function( + p, df, point.shape = 19, point.alpha = 0.4, ...){ + p <- .add_points_layer( + p, df, + x = attributes(df)[["x"]], + y = attributes(df)[["y"]], + point.shape = point.shape, + point.alpha = point.alpha, + ...) + return(p) +} + +# This function adds ellipse visualization. +.add_ellipse_to_ordination <- function( + p, df, + ellipse.alpha = 0.2, + ellipse.linewidth = if(is.null(attributes(df)[["fill.by"]])) 0.5 else 0, + ellipse.linetype = 1, + confidence.level = 0.95, + ...){ + # To disable "no visible binding for global variable" message in cmdcheck + color <- NULL + # + if( !.are_whole_numbers(ellipse.linetype) ){ + stop("'vec.linetype' must be a whole number.", call. = FALSE) + } + if ( !(is.numeric(ellipse.alpha) && ellipse.alpha > 0 && + ellipse.alpha < 1 ) ) { + stop("'ellipse.alpha' must be a number between 0 and 1.", call. = FALSE) + } + if ( !(is.numeric(ellipse.linewidth) && ellipse.linewidth >= 0) ) { + stop("'ellipse.linewidth' must be a positive number.", call. = FALSE) + } + if( !(is.numeric(confidence.level) && confidence.level > 0 && + confidence.level < 1) ) { + stop("'confidence.level' must be a number between 0 and 1.", + call. = FALSE) + } + # + # Get all the arguments. User can fill and colour the ellipses separately. + # However, in most of the cases that might not make sense, but it is still + # made possible. + args <- list( + mapping = aes( + group = if( !is.null(attributes(df)[["group.by"]]) ) + .data[[attributes(df)[["group.by"]]]], + colour = if( !is.null(attributes(df)[["colour.by"]]) ) + .data[[attributes(df)[["colour.by"]]]], + fill = if( !is.null(attributes(df)[["fill.by"]]) ) + .data[[attributes(df)[["fill.by"]]]], + ), + geom = "polygon", + linewidth = ellipse.linewidth, + linetype = ellipse.linetype, + level = confidence.level, + alpha = if( is.null(attributes(df)[["fill.by"]]) ) 0 else ellipse.alpha + ) + # If user did not specify coloring, add black border to ellipses + if( is.null(attributes(df)[["colour.by"]]) ){ + args[["color"]] <- "black" + } + p <- p + do.call(stat_ellipse, args) + return(p) +} + +# This method connects points with a line or arrow. +.connect_ordination_points <- function(p, df){ + # Get arguments. If user sorted the data, use directed arrow. + args <- list( + mapping = aes( + x = .data[[attributes(df)[["x"]]]], + y = .data[[attributes(df)[["y"]]]], + group = .data[[attributes(df)[["pair.by"]]]], + ), + arrow = if(!is.null(attributes(df)[["sort.by"]])) + arrow(length = unit(0.2, "cm"), type = "closed"), + alpha = 0.4 + ) + args <- args[ lengths(args) > 0 ] + p <- p + do.call(geom_path, args) + return(p) +} + +# This method adds group +.add_centroids <- function(p, df, ...){ + df_centroids <- attributes(df)[["centroids"]] + # Add centroids + p <- p + + geom_point( + data = df_centroids, + aes(x = x_centroid, y = y_centroid), + size = 5, + shape = 4, # cross + stroke = 2 # thick border + ) + return(p) +} + +# This method adds group +.add_centroid_lines <- function(p, df, grouping_var, ...){ + df_centroids <- attributes(df)[["centroids"]] + # Add centroids data to original df + df <- df %>% + dplyr::left_join(df_centroids, by = grouping_var) + # Connect points with centroids + p <- p + geom_segment(data = df, aes( + x = x_centroid, y = x_centroid, + xend = .data[[attributes(df)[["x"]]]], + yend = .data[[attributes(df)[["y"]]]] + ), alpha = 0.4) + return(p) +} + +# This methods creates vectors that start from global mean and ends to group +# centroids. This shows how the covariate correlates with the ordination. +.add_centroids_vector <- function(p, df, ...){ + # Calculate centroids + grouping_var <- c( + attributes(df)[["group.by"]], attributes(df)[["fill.by"]]) |> unique() + df_centroids <- df |> + group_by(across(all_of(grouping_var))) |> + summarise( + x_centroid = mean(.data[[attributes(df)[["x"]]]], na.rm = TRUE), + y_centroid = mean(.data[[attributes(df)[["y"]]]], na.rm = TRUE), + .groups = "drop" + ) + # Add global mean + df_centroids[["x_global"]] <- mean( + df[[attributes(df)[["x"]]]], na.rm = TRUE) + df_centroids[["y_global"]] <- mean( + df[[attributes(df)[["y"]]]], na.rm = TRUE) + # Visualize vectors + p <- p + geom_segment( + data = df_centroids, + aes( + x = x_global, y = y_global, + xend = x_centroid, yend = x_centroid, + ), + arrow = arrow(length = unit(0.2, "cm"), type = "closed"), + size = 1 + ) + # Add label to denote which vector belongs to which group + p <- p + ggrepel::geom_label_repel( + data = df_centroids, + mapping = aes( + x = x_centroid, + y = y_centroid, + label = .data[[grouping_var]]) + ) + return(p) +} + +# This methods adds "species scores", i.e., the coordinates of species in the +# ordination. +.add_rotation <- function(p, df){ + # The points are added only if the rotation is present in the data + rot_names <- c("rotation") + if( !is.null(attributes(df)[["rotation"]]) ){ + df_species <- attributes(df)[["rotation"]] + # Add points + p <- p + geom_point( + data = df_species, + mapping = aes( + x = .data[[attributes(df)[["x"]]]], + y = .data[[attributes(df)[["y"]]]] + ), + colour = "red", size = 1) + } + return(p) +} + +# This methods coors the backrounds based on the point density. This creates +# "landscape plots". +.add_point_density <- function(p, df, adjust = 1, ...){ + # Calculate the correct smoothing bandwidth. We could use the default + # values, but this custom bandwidth helps us to highlight the density + # better. + bandwidth <- adjust * c( + .get_bandwidth(df[[attributes(df)[["x"]]]]), + .get_bandwidth(df[[attributes(df)[["y"]]]])) + # Add the landscape + p <- p + + stat_density_2d( + aes(fill = after_stat(density)), + geom = "raster", + h = bandwidth, + contour = FALSE, + alpha = 0.5, + ) + + scale_fill_gradient(name = "Density", low = "white", high = "black") + return(p) +} + +# This function calculates the +.get_bandwidth <- function(x){ + r <- quantile(x, c(0.25, 0.75)) + 4 * 1.06 * min(sd(x), (r[[2]] - r[[1]])/1.34) * length(x)^(-.2) +} + +# This function adjusts the theme of the plot +.adjust_ordination_theme <- function(p, df, ...){ p <- p + theme_classic() p <- p + labs(x = attributes(df)[["x"]]) if( !is.null(attributes(df)[["fill.by"]]) ){ @@ -190,4 +519,4 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), p <- p + labs(shape = attributes(df)[["size.by"]]) } return(p) -} +} \ No newline at end of file From 214f24b819c8e6c5842c679329b0547838b66bf7 Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Mon, 12 May 2025 16:46:29 +0300 Subject: [PATCH 06/28] up --- R/plotOrdination.R | 73 ++++++++++++++++++++++++++++++---------------- 1 file changed, 48 insertions(+), 25 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index ce3b3e34..8e7af440 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -63,7 +63,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), assay.type = "counts", add.points = TRUE, add.ellipse = FALSE, add.density = FALSE, add.centroids = FALSE, add.centroids.lines = FALSE, add.vectors = FALSE, - add.rotation = FALSE, + add.rotation = FALSE, add.expl.var = FALSE, ...){ # Check if there are any reduced dim present if( length(reducedDims(x)) == 0L ){ @@ -125,6 +125,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !.is_a_bool(add.rotation) ){ stop("'add.rotation' must be TRUE or FALSE.", call. = FALSE) } + if( !.is_a_bool(add.expl.var) ){ + stop("'add.expl.var' must be TRUE or FALSE.", call. = FALSE) + } # add.density cannot be specified simultaneously with fill.by as they both # are using fill aesthetic and we can have only on fill scale. @@ -144,7 +147,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), pair.by = pair.by, sort.by = sort.by, facet.by = facet.by, - assay.type = "counts" + assay.type = assay.type, + add.expl.var = add.expl.var ) args <- c(args, list(...)) return(args) @@ -153,7 +157,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # This function retrieves the data from reducedDim .get_ordination_data <- function( x, dimred, ncomponents, colour.by, fill.by, shape.by, size.by, group.by, - pair.by, sort.by, facet.by, assay.type, ...){ + pair.by, sort.by, facet.by, assay.type, add.expl.var = FALSE, ...){ # Get data and store the original attributes that might include rotation # data, for instance df <- reducedDim(x, dimred) @@ -181,6 +185,26 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), colnames(rotation) <- colnames(df) } + expl_var_name <- c("eig") + xlab <- x_var + ylab <- y_var + if( add.expl.var && any(expl_var_name %in% names(orig_attributes)) ){ + eigen <- orig_attributes[expl_var_name][[1L]] + xlab <- paste0( + xlab, " (", + round(eigen[ncomponents][[1L]], 1), + "%)" + ) + ylab <- paste0( + ylab, " (", + round(eigen[ncomponents][[2L]], 1), + "%)" + ) + } else if( add.expl.var ){ + warning("No explained variance found from the data.", call. = FALSE) + } + + # List data that is fetched from colData cols <- c(shape.by, size.by, group.by, pair.by, sort.by, facet.by) # colour.by and fill.by can also specify a feature and its abundance. @@ -234,7 +258,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), pair.by = pair.by, sort.by = sort.by, facet.by = facet.by, - assay.type = assay.type + assay.type = assay.type, + xlab = xlab, + ylab = ylab ) attr(df, "rotation") <- rotation attr(df, "centroids") <- df_centroids @@ -287,7 +313,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Add vectors from global centroid to group centroids if( add.vectors && length(grouping_var) > 0L ){ - p <- .add_centroids_vector(p, df, ...) + p <- .add_centroids_vector(p, df, grouping_var, ...) } # If facetting was specified, split plot to separate panels if( !is.null(attributes(df)[["facet.by"]]) ){ @@ -409,7 +435,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), dplyr::left_join(df_centroids, by = grouping_var) # Connect points with centroids p <- p + geom_segment(data = df, aes( - x = x_centroid, y = x_centroid, + x = x_centroid, y = y_centroid, xend = .data[[attributes(df)[["x"]]]], yend = .data[[attributes(df)[["y"]]]] ), alpha = 0.4) @@ -418,28 +444,14 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # This methods creates vectors that start from global mean and ends to group # centroids. This shows how the covariate correlates with the ordination. -.add_centroids_vector <- function(p, df, ...){ - # Calculate centroids - grouping_var <- c( - attributes(df)[["group.by"]], attributes(df)[["fill.by"]]) |> unique() - df_centroids <- df |> - group_by(across(all_of(grouping_var))) |> - summarise( - x_centroid = mean(.data[[attributes(df)[["x"]]]], na.rm = TRUE), - y_centroid = mean(.data[[attributes(df)[["y"]]]], na.rm = TRUE), - .groups = "drop" - ) - # Add global mean - df_centroids[["x_global"]] <- mean( - df[[attributes(df)[["x"]]]], na.rm = TRUE) - df_centroids[["y_global"]] <- mean( - df[[attributes(df)[["y"]]]], na.rm = TRUE) +.add_centroids_vector <- function(p, df, grouping_var, ...){ + df_centroids <- attributes(df)[["centroids"]] # Visualize vectors p <- p + geom_segment( data = df_centroids, aes( x = x_global, y = y_global, - xend = x_centroid, yend = x_centroid, + xend = x_centroid, yend = y_centroid, ), arrow = arrow(length = unit(0.2, "cm"), type = "closed"), size = 1 @@ -503,9 +515,20 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # This function adjusts the theme of the plot -.adjust_ordination_theme <- function(p, df, ...){ +.adjust_ordination_theme <- function( + p, df, + xlab = attributes(df)[["xlab"]], + ylab = attributes(df)[["ylab"]], + ...){ + if( !.is_a_string(xlab) ){ + stop("'xlab' must be a single character value.", call. = FALSE) + } + if( !.is_a_string(ylab) ){ + stop("'ylab' must be a single character value.", call. = FALSE) + } + # p <- p + theme_classic() - p <- p + labs(x = attributes(df)[["x"]]) + p <- p + labs(x = xlab, y = ylab) if( !is.null(attributes(df)[["fill.by"]]) ){ p <- p + labs(fill = attributes(df)[["fill.by"]]) } From ff1b28d71004787a8ade9548eb9b2e29abf2dccd Mon Sep 17 00:00:00 2001 From: TuomasBorman Date: Thu, 15 May 2025 09:08:50 +0300 Subject: [PATCH 07/28] up --- R/plotOrdination.R | 44 ++++++++++++++++++++++++++++++++++++-------- 1 file changed, 36 insertions(+), 8 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 8e7af440..90b74458 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -57,6 +57,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), shape.by = NULL, size.by = NULL, group.by = NULL, + linetype.by = NULL, pair.by = NULL, sort.by = NULL, facet.by = NULL, @@ -94,6 +95,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), temp <- .check_metadata_variable(tse, fill.by, FALSE, TRUE, FALSE, TRUE) temp <- .check_metadata_variable(tse, shape.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, size.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable( + tse, linetype.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, group.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, pair.by, FALSE, TRUE, FALSE, FALSE) temp <- .check_metadata_variable(tse, sort.by, FALSE, TRUE, FALSE, FALSE) @@ -144,6 +147,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), shape.by = shape.by, size.by = size.by, group.by = group.by, + linetype.by = linetype.by, pair.by = pair.by, sort.by = sort.by, facet.by = facet.by, @@ -157,7 +161,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # This function retrieves the data from reducedDim .get_ordination_data <- function( x, dimred, ncomponents, colour.by, fill.by, shape.by, size.by, group.by, - pair.by, sort.by, facet.by, assay.type, add.expl.var = FALSE, ...){ + linetype.by, pair.by, sort.by, facet.by, assay.type, + add.expl.var = FALSE, ...){ # Get data and store the original attributes that might include rotation # data, for instance df <- reducedDim(x, dimred) @@ -204,9 +209,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), warning("No explained variance found from the data.", call. = FALSE) } - # List data that is fetched from colData - cols <- c(shape.by, size.by, group.by, pair.by, sort.by, facet.by) + cols <- c( + shape.by, size.by, group.by, linetype.by, pair.by, sort.by, facet.by) # colour.by and fill.by can also specify a feature and its abundance. # Either get their column name to fetch from colData or directly add them # to data. @@ -225,6 +230,20 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), cd <- colData(x)[, cols, drop = FALSE] df <- cbind(df, cd) + # Check that the values are correct + if( !is.null(group.by) && is.numeric(df[[group.by]]) ){ + stop("Values specified by 'group.by' must be categorical.", + call. = FALSE) + } + if( !is.null(linetype.by) && is.numeric(df[[linetype.by]]) ){ + stop("Values specified by 'linetype.by' must be categorical.", + call. = FALSE) + } + if( !is.null(facet.by) && is.numeric(df[[facet.by]]) ){ + stop("Values specified by 'facet.by' must be categorical.", + call. = FALSE) + } + # Sort the data. For instance, if we want to add arrows between consecutive # time points, this is essential step. if( !is.null(sort.by) ){ @@ -255,6 +274,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), shape.by = shape.by, size.by = size.by, group.by = group.by, + linetype.by = linetype.by, pair.by = pair.by, sort.by = sort.by, facet.by = facet.by, @@ -379,16 +399,20 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), .data[[attributes(df)[["colour.by"]]]], fill = if( !is.null(attributes(df)[["fill.by"]]) ) .data[[attributes(df)[["fill.by"]]]], + linetype = if( !is.null(attributes(df)[["linetype.by"]]) ) + .data[[attributes(df)[["linetype.by"]]]], ), geom = "polygon", linewidth = ellipse.linewidth, - linetype = ellipse.linetype, level = confidence.level, alpha = if( is.null(attributes(df)[["fill.by"]]) ) 0 else ellipse.alpha ) # If user did not specify coloring, add black border to ellipses if( is.null(attributes(df)[["colour.by"]]) ){ - args[["color"]] <- "black" + + } + if( is.null(attributes(df)[["linetype.by"]]) ){ + args[["linetype"]] = ellipse.linetype } p <- p + do.call(stat_ellipse, args) return(p) @@ -508,10 +532,11 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), return(p) } -# This function calculates the +# This function calculates the smoothing bandwidth. It highlights better point- +# rich areas than the default choice. .get_bandwidth <- function(x){ r <- quantile(x, c(0.25, 0.75)) - 4 * 1.06 * min(sd(x), (r[[2]] - r[[1]])/1.34) * length(x)^(-.2) + 4 * 1.06 * min(sd(x), (r[[2]] - r[[1]])/1.34) * length(x)^(-0.2) } # This function adjusts the theme of the plot @@ -541,5 +566,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !is.null(attributes(df)[["size.by"]]) ){ p <- p + labs(shape = attributes(df)[["size.by"]]) } + if( !is.null(attributes(df)[["linetype.by"]]) ){ + p <- p + labs(linetype = attributes(df)[["linetype.by"]]) + } return(p) -} \ No newline at end of file +} From 3435451bd0b7425695af2df86994f4b42cf51486 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Sat, 28 Feb 2026 21:01:32 +0200 Subject: [PATCH 08/28] Equal scale for axes --- R/plotOrdination.R | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 90b74458..7d914b50 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -345,6 +345,11 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Adjust theme p <- .adjust_ordination_theme(p, df, ...) + + # Enforce same scale to x and y axis. Without equal scale, the results and + # interpretations might be misleading + p <- p + coord_equal() + return(p) } From d73aedac653f13b5dee1bac82628774a1e3232e6 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 22:45:24 +0300 Subject: [PATCH 09/28] up --- DESCRIPTION | 4 +- NAMESPACE | 8 +- R/AllGenerics.R | 5 + R/plotCCA.R | 267 ++------------------- R/plotJointRPCA.R | 340 +++++++++++++++++++++++++++ R/plotOrdination.R | 304 ++++++++++++++++++++---- R/utils.R | 2 + man/plotCCA.Rd | 35 +-- man/plotJointRPCA.Rd | 121 ++++++++++ man/plotOrdination.Rd | 253 ++++++++++++++++++++ tests/testthat/test-plotJointRPCA.R | 287 ++++++++++++++++++++++ tests/testthat/test-plotOrdination.R | 228 ++++++++++++++++++ 12 files changed, 1536 insertions(+), 318 deletions(-) create mode 100644 R/plotJointRPCA.R create mode 100644 man/plotJointRPCA.Rd create mode 100644 man/plotOrdination.Rd create mode 100644 tests/testthat/test-plotJointRPCA.R create mode 100644 tests/testthat/test-plotOrdination.R diff --git a/DESCRIPTION b/DESCRIPTION index 117fadfa..a991d385 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -48,6 +48,7 @@ Imports: ggnewscale, ggrepel, ggtree, + grid, methods, rlang, S4Vectors, @@ -77,8 +78,7 @@ Suggests: Remotes: github::microbiome/miaTime Roxygen: list(markdown = TRUE) -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 VignetteBuilder: knitr URL: https://github.com/microbiome/miaViz BugReports: https://github.com/microbiome/miaViz/issues - diff --git a/NAMESPACE b/NAMESPACE index aee57831..cef43b04 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -9,9 +9,11 @@ export(plotColTile) export(plotDMNFit) export(plotFeaturePrevalence) export(plotHistogram) +export(plotJointRPCA) export(plotLoadings) export(plotMediation) export(plotNMDS) +export(plotOrdination) export(plotPrevalence) export(plotPrevalentAbundance) export(plotRDA) @@ -37,8 +39,10 @@ exportMethods(plotColTree) exportMethods(plotDMNFit) exportMethods(plotFeaturePrevalence) exportMethods(plotHistogram) +exportMethods(plotJointRPCA) exportMethods(plotLoadings) exportMethods(plotMediation) +exportMethods(plotOrdination) exportMethods(plotPrevalence) exportMethods(plotPrevalentAbundance) exportMethods(plotRDA) @@ -100,6 +104,7 @@ importFrom(dplyr,rename) importFrom(dplyr,rename_with) importFrom(dplyr,row_number) importFrom(dplyr,select) +importFrom(dplyr,slice_max) importFrom(dplyr,summarise) importFrom(dplyr,summarize) importFrom(dplyr,ungroup) @@ -146,13 +151,14 @@ importFrom(ggtree,ggtree) importFrom(ggtree,gheatmap) importFrom(ggtree,groupOTU) importFrom(ggtree,theme_tree) +importFrom(grid,arrow) +importFrom(grid,unit) importFrom(mia,meltSE) importFrom(rlang,"!!") importFrom(rlang,":=") importFrom(rlang,sym) importFrom(scales,pretty_breaks) importFrom(scales,rescale) -importFrom(scater,plotReducedDim) importFrom(scater,retrieveCellInfo) importFrom(scater,retrieveFeatureInfo) importFrom(stats,formula) diff --git a/R/AllGenerics.R b/R/AllGenerics.R index 87d2e071..e5431307 100644 --- a/R/AllGenerics.R +++ b/R/AllGenerics.R @@ -155,3 +155,8 @@ setGeneric("plotBoxplot", signature = c("object"), function(object, ...) #' @export setGeneric("plotOrdination", signature = c("x"), function(x, ...) standardGeneric("plotOrdination")) + +#' @rdname plotJointRPCA +#' @export +setGeneric("plotJointRPCA", signature = c("x"), function(x, ...) + standardGeneric("plotJointRPCA")) diff --git a/R/plotCCA.R b/R/plotCCA.R index 9c5a7494..b903f912 100644 --- a/R/plotCCA.R +++ b/R/plotCCA.R @@ -6,8 +6,7 @@ #' #' @param x a #' \code{\link[TreeSummarizedExperiment:TreeSummarizedExperiment-constructor]{TreeSummarizedExperiment}} -#' or a matrix of weights. The latter is returned as output from -#' \code{\link[mia:runCCA]{getRDA}}. +#' object. #' #' @param dimred \code{Character scalar} or \code{integer scalar}. Determines #' the reduced dimension to @@ -120,14 +119,17 @@ #' formula = assay ~ ClinicalStatus + Gender + Age, #' distance = "bray", #' na.action = na.exclude -#' ) +#' ) #' #' suppressWarnings({ #' # Create RDA plot coloured by variable #' plotRDA(tse, "RDA", colour.by = "ClinicalStatus") #' -#' # Create RDA plot with empty ellipses -#' plotRDA(tse, "RDA", colour.by = "ClinicalStatus", add.ellipse = "colour") +#' # Create RDA plot with ellipses +#' plotRDA( +#' tse, "RDA", colour.by = "ClinicalStatus", fill.by = "ClinicalStatus", +#' add.ellipse = TRUE +#' ) #' #' # Create RDA plot with text encased in labels #' plotRDA(tse, "RDA", colour.by = "ClinicalStatus", vec.text = FALSE) @@ -137,20 +139,13 @@ #' #' # Create RDA plot without vectors #' plotRDA(tse, "RDA", colour.by = "ClinicalStatus", add.vectors = FALSE) -#' -#' # Calculate RDA as a separate object -#' rda_mat <- getRDA( -#' tse, -#' assay.type = "relabundance", -#' formula = assay ~ ClinicalStatus + Gender + Age, -#' distance = "bray", -#' na.action = na.exclude -#' ) -#' -#' # Create RDA plot from RDA matrix -#' plotRDA(rda_mat) #' }) #' +#' @seealso +#' \itemize{ +#' \item \code{\link[=plotOrdination]{plotOrdination}} +#' } +#' NULL #' @rdname plotCCA @@ -163,16 +158,6 @@ setMethod("plotCCA", signature = c(x = "SingleCellExperiment"), } ) -#' @rdname plotCCA -#' @aliases plotRDA -#' @export -setMethod("plotCCA", signature = c(x = "matrix"), - function(x, ...){ - # Reproduce plotRDA function - return(plotRDA(x, ...)) - } -) - #' @rdname plotCCA #' @aliases plotCCA #' @export @@ -201,26 +186,10 @@ setMethod("plotRDA", signature = c(x = "SingleCellExperiment"), args[["ncomponents"]] <- 2L # Get data for plotting plot_args <- list() - plot_args[["ellipse_data"]] <- do.call(.get_rda_ellipse_data, args) plot_args[["vector_data"]] <- do.call(.get_rda_vector_data, args) - plot_args[["centroids"]] <- do.call(.get_rda_centroids_data, args) - plot_args[["species_scores"]] <- do.call(.get_rda_species_data, args) - plot_args[["plot"]] <- do.call(.create_rda_baseplot, args) + p <- plotOrdination(x, dimred = dimred, ...) # Create a final plot - p <- .rda_plotter(plot_args, ...) - return(p) - } -) - -#' @rdname plotCCA -#' @aliases plotCCA -#' @export -setMethod("plotRDA", signature = c(x = "matrix"), - function(x, ...){ - # Construct TreeSE from rda/cca object - x <- .rda2tse(x) - # Run plotRDA method for TreeSE - p <- plotRDA(x, "RDA", ...) + p <- .rda_plotter_vector(p, plot_args, ...) return(p) } ) @@ -269,34 +238,6 @@ setMethod("plotRDA", signature = c(x = "matrix"), return(reduced_dim) } -# This function retrieves optional data that is used for creating an ellipses. -#' @importFrom scater retrieveCellInfo -.get_rda_ellipse_data <- function( - tse, reduced_dim, add.ellipse = TRUE, colour_by = color_by, - color_by = colour.by, colour.by = color.by, color.by = NULL, ...){ - # - if( !(add.ellipse %in% c(TRUE, FALSE, "fill", "color", "colour") && - length(add.ellipse) == 1L ) ){ - stop("'add.ellipse' must be one of c(TRUE, FALSE, 'fill', ", - "'color').", call. = FALSE) - } - if( !(is.null(colour_by) || .is_a_string(colour_by) && - colour_by %in% colnames(colData(tse)) ) ){ - stop("'colour_by' must be NULL or name of column from colData(x).", - call. = FALSE) - } - # - ellipse_data <- NULL - if( add.ellipse != FALSE && !is.null(colour_by) ){ - # Ellipse data is the same ordination data - ellipse_data <- as.data.frame(reduced_dim) - # Add sample metadata from colData - ellipse_data[[colour_by]] <- retrieveCellInfo(tse, colour_by)[["value"]] - attributes(ellipse_data)[["colour_by"]] <- colour_by - } - return(ellipse_data) -} - # This function retrieves data for creating vectors. Moreover, it wrangles the # vector data and controls what information is added to vector text or labels. .get_rda_vector_data <- function( @@ -487,43 +428,6 @@ setMethod("plotRDA", signature = c(x = "matrix"), return(vector_data) } -# This functions returns optional centroids for plotting. -.get_rda_centroids_data <- function( - reduced_dim, add.centroids = FALSE, ncomponents = 2L, ...){ - # - if( !.is_a_bool(add.centroids) ){ - stop("'add.centroids' must be TRUE or FALSE.", call. = FALSE) - } - if( !.is_an_integer(ncomponents) ){ - stop("'ncomponents' must be an integer.", call. = FALSE) - } - # - res <- if(add.centroids) .get_rda_attribute(reduced_dim, "centroids") - if( !is.null(res) ){ - res <- res[, seq_len(ncomponents), drop = FALSE] |> as.data.frame() - colnames(res) <- c("x", "y") - } - return(res) -} - -# This functions returns optional species scores for plotting. -.get_rda_species_data <- function( - reduced_dim, add.species = FALSE, ncomponents = 2L, ...){ - if( !.is_a_bool(add.species) ){ - stop("'add.species' must be TRUE or FALSE.", call. = FALSE) - } - if( !.is_an_integer(ncomponents) ){ - stop("'ncomponents' must be an integer.", call. = FALSE) - } - # - res <- if(add.species) .get_rda_attribute(reduced_dim, "species") - if( add.species ){ - res <- res[, seq_len(ncomponents), drop = FALSE] |> as.data.frame() - colnames(res) <- c("x", "y") - } - return(res) -} - # This function is used to fetch specified datatype from attributes if it # exists. .get_rda_attribute <- function(reduced_dim, attr_names){ @@ -536,127 +440,6 @@ setMethod("plotRDA", signature = c(x = "matrix"), return(res) } -# This function utilizes scater::plotReducedDim to create "baseplot". To where -# we can build the the plot. The idea is that the theme is similar in all -# ordination plots. -#' @importFrom scater plotReducedDim -.create_rda_baseplot <- function( - tse, dimred, reduced_dim, ncomponents = 2L, - add.expl.var = FALSE, expl.var = expl_var, expl_var = NULL, - colour_by = color_by, color_by = colour.by, - colour.by = color.by, color.by = NULL, ...){ - # - if( !.is_a_bool(add.expl.var) ){ - stop("'add.expl.var' must be TRUE or FALSE.", call. = FALSE) - } - if( !.is_an_integer(ncomponents) ){ - stop("'ncomponents' must be an integer.", call. = FALSE) - } - if( !( is.null(expl.var) || (is.numeric(expl.var) && - length(expl.var) == ncomponents )) ){ - stop("'expl.var' must be numeric vector with length ", ncomponents, - ".", call. = FALSE) - } - if( !(is.null(colour_by) || .is_a_string(colour_by) && - colour_by %in% colnames(colData(tse)) ) ){ - stop("'colour_by' must be NULL or name of column from colData(x).", - call. = FALSE) - } - # - # If specified, get explained variance - if( add.expl.var && is.null(expl.var) ){ - eigen_vals <- attr(reduced_dim, "eig") - # Convert to explained variance and take only first two components - expl_var <- eigen_vals / sum(eigen_vals) - expl_var <- expl_var[seq_len(ncomponents)]*100 - } - # Create argument list - args <- c(list(object = tse, dimred = dimred, ncomponents = ncomponents, - colour_by = colour_by, percentVar = expl_var), list(...)) - # Remove additional arguments since plotReducedDim fails if we feed - # values that are not recognized - remove <- names(args) %in% c( - "add.significance", "add.expl.var", "add.ellipse", "add.vectors", - "vec.lab", "sep.group", "repl.underscore", "add.centroids", - "add.species", "ellipse.alpha", "ellipse.linewidth", "ellipse.linetype", - "confidence.level", "vec.size", "vec.color", "vec.colour", - "vec.linetype", "arrow.size", "min.segment.length", "label.color", - "label.colour", "label.size", "parse.labels", "vec.text", - "repel.labels", "position", "nudge_x", "nudge_y", "direction", - "max.overlaps", "check_overlap", "ignore.case") - args <- args[ !remove ] - # Get scatter plot with plotReducedDim --> keep theme similar between - # ordination methods - p <- do.call(plotReducedDim, args) - return(p) -} - -# This function is used to create the plot. -.rda_plotter <- function(plot_data, ...){ - # Get the scatter plot - plot <- plot_data[["plot"]] - # Add ellipse - plot <- .rda_plotter_ellipse(plot, plot_data, ...) - # Add vectors - plot <- .rda_plotter_vector(plot, plot_data, ...) - # Add centroids - plot <- .rda_plotter_centroids_or_species(plot, plot_data, "centroids") - # Add species - plot <- .rda_plotter_centroids_or_species(plot, plot_data, "species_scores") - return(plot) -} - -# This function adds ellipse visualization. -.rda_plotter_ellipse <- function( - plot, plot_data, add.ellipse = TRUE, ellipse.alpha = 0.2, - ellipse.linewidth = 0.1, ellipse.linetype = 1, confidence.level = 0.95, - ...){ - # To disable "no visible binding for global variable" message in cmdcheck - color <- NULL - # - if( !(add.ellipse %in% c(TRUE, FALSE, "fill", "color", "colour") && - length(add.ellipse) == 1L ) ){ - stop("'add.ellipse' must be one of c(TRUE, FALSE, 'fill', ", - "'color').", call. = FALSE) - } - if( !.are_whole_numbers(ellipse.linetype) ){ - stop("'vec.linetype' must be a whole number.", call. = FALSE) - } - if ( !(is.numeric(ellipse.alpha) && ellipse.alpha > 0 && - ellipse.alpha < 1 ) ) { - stop("'ellipse.alpha' must be a number between 0 and 1.", call. = FALSE) - } - if ( !(is.numeric(ellipse.linewidth) && ellipse.linewidth > 0) ) { - stop("'ellipse.linewidth' must be a positive number.", call. = FALSE) - } - if( !(is.numeric(confidence.level) && confidence.level > 0 && - confidence.level < 1) ) { - stop("'confidence.level' must be a number between 0 and 1.", - call. = FALSE) - } - # - data <- plot_data[["ellipse_data"]] - if( !is.null(data) ){ - xvar <- colnames(data)[[1]] - yvar <- colnames(data)[[2]] - colour_var <- attributes(data)[["colour_by"]] - # Add ellipses to plot (fill or colour the edge) - fill <- add.ellipse %in% c(TRUE, "fill") - plot <- plot + stat_ellipse( - data = data, - mapping = aes( - x = .data[[xvar]], y = .data[[yvar]], - color = .data[[colour_var]], fill = after_scale(color)), - geom = "polygon", - linewidth = ellipse.linewidth, - linetype = ellipse.linetype, - level = confidence.level, - alpha = if(fill) ellipse.alpha else 0 - ) - } - return(plot) -} - # This function adds vector and text layer to the plot. #' @importFrom ggrepel geom_text_repel geom_label_repel .rda_plotter_vector <- function( @@ -666,7 +449,7 @@ setMethod("plotRDA", signature = c(x = "matrix"), vec.text = TRUE, repel.labels = TRUE, parse.labels = TRUE, add.significance = TRUE, vec.linetype = 1, min.segment.length = 5, position = NULL, nudge_x = NULL, nudge_y = NULL, direction = "both", - max.overlaps = 10, check_overlap = FALSE, ... + max.overlaps = 10, ... ){ # if ( !(is.numeric(vec.size) && vec.size > 0) ) { @@ -714,7 +497,7 @@ setMethod("plotRDA", signature = c(x = "matrix"), x = 0, y = 0, xend = .data[[xvar]], yend = .data[[yvar]], group = .data[["group"]]), arrow = arrow(length = unit(arrow.size, "cm")), - color = vec.color, linetype = vec.linetype, size = vec.size) + color = vec.color, linetype = vec.linetype, linewidth = vec.size) # Add vector labels (text or label) # Make list of arguments for geom_text/geom_label label_args <- list( @@ -734,8 +517,6 @@ setMethod("plotRDA", signature = c(x = "matrix"), max.time = 0.5, max.iter = 10000, max.overlaps = max.overlaps, direction = direction, seed = NA, verbose = FALSE ) - } else if( !repel.labels && vec.text ){ - label_args <- c(label_args, check_overlap = check_overlap) } # Choose right function and call it FUN <- if( repel.labels && vec.text ) geom_text_repel @@ -746,19 +527,3 @@ setMethod("plotRDA", signature = c(x = "matrix"), } return(plot) } - -# This function adds centroids or species layer to the plot. -.rda_plotter_centroids_or_species <- function(plot, plot_data, type){ - # To disable "no visible binding for global variable" message in cmdcheck - x <- y <- NULL - data <- plot_data[[type]] - if( !is.null(data) ){ - plot <- plot + geom_point( - data, - mapping = aes(x = x, y = y), - shape = if(type == "centroids") 10L else 4L, - color = if(type == "centroids") "blue" else "red", - ) - } - return(plot) -} diff --git a/R/plotJointRPCA.R b/R/plotJointRPCA.R new file mode 100644 index 00000000..38aa222b --- /dev/null +++ b/R/plotJointRPCA.R @@ -0,0 +1,340 @@ +#' @name +#' plotJointRPCA +#' +#' @title +#' Visualize Joint-RPCA results +#' +#' @description +#' Creates a two-dimensional ordination plot from Joint-RPCA results stored in a +#' \code{SingleCellExperiment} or \code{MultiAssayExperiment}. In addition to the +#' sample ordination, feature loadings can be visualized as vectors. +#' +#' @details +#' This function is a wrapper around \code{\link[=plotOrdination]{plotOrdination}} +#' for Joint-RPCA results. Consequently, most graphical parameters are passed +#' directly to \code{plotOrdination()}, including options for colouring, +#' grouping, faceting and adding ellipses. See +#' \code{\link[=plotOrdination]{plotOrdination}} for a complete description of +#' these arguments. +#' +#' Feature loadings are plotted as vectors originating from the origin. By +#' default, only the longest loading vectors are shown for each data layer to +#' improve readability. +#' +#' @return +#' A \code{ggplot2} object. +#' +#' @param x A +#' \code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}} +#' or +#' \code{\link[MultiAssayExperiment:MultiAssayExperiment-class]{MultiAssayExperiment}} +#' object containing Joint-RPCA results. +#' +#' @param dimred \code{character(1)}: Specifies the name of the Joint-RPCA +#' result to plot. +#' +#' @param add.vectors: \code{logical(1)} or \code{character}. If +#' \code{TRUE}, feature loading vectors are added. Alternatively, a character +#' vector can be supplied to display only features whose names match the given +#' pattern(s). (Default: \code{TRUE}) +#' +#' @param ... +#' Additional arguments passed to +#' \code{\link[=plotOrdination]{plotOrdination}} and to the feature-vector +#' plotting. +#' Commonly used Joint-RPCA-specific arguments include: +#' \itemize{ +#' \item \code{ntop}: \code{integer(1)} or \code{NULL}. Maximum number of +#' loading vectors shown per data layer. The vectors with the largest lengths +#' are retained. Set to \code{NULL} to display all vectors. (Default: \code{10}) +#' } +#' +#' @examples +#' data("HintikkaXOData") +#' mae <- HintikkaXOData +#' +#' mae[[1]] <- transformAssay( +#' mae[[1]], +#' assay.type = "counts", +#' method = "rclr", +#' impute = FALSE +#' ) +#' mae[[2]] <- transformAssay( +#' mae[[2]], +#' assay.type = "nmr", +#' method = "log10" +#' ) +#' +#' mae <- addJointRPCA( +#' mae, +#' experiments = c(1, 2), +#' assay.types = c("rclr", "log10") +#' ) +#' +#' # Basic Joint-RPCA plot +#' plotJointRPCA(mae, "JointRPCA") +#' +#' # Colour samples by metadata +#' plotJointRPCA(mae, "JointRPCA", colour.by = "Fat") +#' +#' # Add confidence ellipses +#' plotJointRPCA( +#' mae, +#' "JointRPCA", +#' colour.by = "Fat", +#' add.ellipse = TRUE +#' ) +#' +#' # Show all loading vectors +#' plotJointRPCA( +#' mae, +#' "JointRPCA", +#' ntop = NULL +#' ) +#' +#' # Display only selected feature vectors +#' plotJointRPCA( +#' mae, +#' "JointRPCA", +#' add.vectors = c("Bacteroides", "Roseburia") +#' ) +#' +#' # Use boxed labels without repelling +#' plotJointRPCA( +#' mae, +#' "JointRPCA", +#' text.labels = FALSE, +#' repel.labels = FALSE +#' ) +#' +#' @seealso +#' \itemize{ +#' \item \code{\link[=plotOrdination]{plotOrdination}} +#' } +#' +NULL + +#' @rdname plotJointRPCA +#' @export +setMethod("plotJointRPCA", signature = c(x = "MultiAssayExperiment"), + function(x, dimred, ...){ + mia:::.check_metadata_present(dimred, x) + # Construct TreeSE from results so that we can use TreeSE function + # (and thus plotOrdination) + res <- metadata(x)[[dimred]] + # + col_data <- x |> colData() + col_data <- col_data[ + match(rownames(res), rownames(col_data)), , drop = FALSE] + rownames(col_data) <- rownames(res) + # + tse <- TreeSummarizedExperiment( + assays = SimpleList(counts = matrix( + ncol = nrow(res), dimnames = list(NULL, rownames(res)))), + colData = col_data, + reducedDims = list(JointRPCA = res) + ) + # Plot + p <- plotJointRPCA(tse, dimred = "JointRPCA", ...) + return(p) + } +) + +#' @rdname plotJointRPCA +#' @export +setMethod("plotJointRPCA", signature = c(x = "SingleCellExperiment"), + function(x, dimred, add.vectors = TRUE, ...){ + .check_dimred_present(dimred, x) + if( !inherits(reducedDim(x, dimred), "JointRPCA") ){ + stop("The plotted results must be a class 'JointRPCA'.", + call. = FALSE) + } + if( !( .is_a_bool(add.vectors) || is.character(add.vectors)) ){ + stop("'add.vectors must be TRUE or FALSE or character vector.", + call. = FALSE) + } + p <- plotOrdination(x, dimred, ...) + vector_data <- .get_joint_rpca_vector_data(x, dimred, add.vectors, ...) + if( !is.null(vector_data ) ){ + p <- .add_ordination_vectors(p, vector_data, ...) + } + return(p) + } +) + +################################ HELP FUNCTIONS ################################ + +# Get feature loadings that will be plotted as vectors +#' @importFrom dplyr group_by slice_max ungroup +.get_joint_rpca_vector_data <- function( + tse, reduced_dim, add.vectors, ignore.case = FALSE, ntop = 10, + ...){ + if( !.is_a_bool(ignore.case) ){ + stop("'ignore.case' must be TRUE or FALSE.", call. = FALSE) + } + if( !(.is_an_integer(ntop) || is.null(ntop)) ){ + stop("'ntop' must be an integer.", call. = FALSE) + } + # + # Get vector data, i.e, rotations or loadings + res <- reducedDim(tse, reduced_dim) + vector_data <- if(!(.is_a_bool(add.vectors) && !add.vectors)) + attributes(res)[["rotation"]] + # If user wanted to plot them, wrangle the data + if( !is.null(vector_data) ){ + # Add labels and layer name + vector_data <- as.data.frame(vector_data) + vector_data[["vector_label"]] <- rownames(vector_data) + n_features <- attr(res, "n_features") + vector_data[["Layer"]] <- rep( + names(n_features), + times = unname(n_features) + ) + + # Subset vectors by selecting only those ones that user has specified + if( is.character(add.vectors) ){ + add.vectors <- paste0(add.vectors, collapse = "|") + keep <- vapply(rownames(vector_data), function(x) + grepl(add.vectors, x, perl = TRUE, ignore.case = ignore.case), + logical(1L)) + vector_data <- vector_data[keep, ] + } + + # Keep only the n longest vectors per layer. The idea is to show only + # features that are highly associated with the showed ordination space. + if (!is.null(ntop) && nrow(vector_data) > 0L) { + xvar <- colnames(vector_data)[1] + yvar <- colnames(vector_data)[2] + + vector_data[["length_temporary"]] <- + sqrt(vector_data[[xvar]]^2 + vector_data[[yvar]]^2) + + vector_data <- vector_data |> + group_by(Layer) |> + slice_max(length_temporary, n = ntop, with_ties = FALSE) |> + ungroup() + + vector_data[["length_temporary"]] <- NULL + } + + # If all vectors were removed, give NULL + if( nrow(vector_data) == 0L ){ + vector_data <- NULL + } + } + + return(vector_data) +} + +# Add feature loadings as vectors +#' @importFrom ggrepel geom_text_repel geom_label_repel +#' @importFrom grid arrow unit +.add_ordination_vectors <- function( + p, + vector_data, + vec.size = 0.5, + arrow.size = 0.25, + label.size = 4, + vec.color = "black", + label.color = "black", + vec.text = TRUE, + text.labels = TRUE, + repel.labels = TRUE, + min.segment.length = 0.5, + box.padding = 0.25, + point.padding = 1e-06, + force = 1, + force_pull = 1, + max.time = 0.5, + max.iter = 10000, + max.overlaps = 10, + direction = "both", + seed = NA, + ... +){ + if ( !(is.numeric(vec.size) && vec.size > 0) ) { + stop("'vec.size' must be a positive number.", call. = FALSE) + } + if ( !(is.numeric(arrow.size) && arrow.size > 0) ) { + stop("'arrow.size' must be a positive number.", call. = FALSE) + } + if ( !(is.numeric(label.size) && label.size > 0) ) { + stop("'label.size' must be a positive number.", call. = FALSE) + } + if ( !.is_non_empty_string(vec.color) ) { + stop("'vec.color' must be a non-empty string specifying a colour", + call. = FALSE) + } + if ( !.is_non_empty_string(label.color) ) { + stop("'label.color' must be a non-empty string specifying a colour", + call. = FALSE) + } + if( !.is_a_bool(vec.text) ){ + stop("'vec.text' must be TRUE or FALSE.", call. = FALSE) + } + if( !.is_a_bool(repel.labels) ){ + stop("'repel.labels' must be TRUE or FALSE.", call. = FALSE) + } + # + # Names of the ordination axes + xvar <- colnames(vector_data)[1] + yvar <- colnames(vector_data)[2] + + # Draw vectors from the origin to the feature coordinates + p <- p + + geom_segment( + data = vector_data, + aes( + x = 0, + y = 0, + xend = .data[[xvar]], + yend = .data[[yvar]], + linetype = .data[["Layer"]] + ), + arrow = arrow(length = unit(arrow.size, "cm")), + linewidth = vec.size, + colour = vec.color, + inherit.aes = FALSE + ) + + # Select the appropriate text/label geometry + FUN <- if (repel.labels) { + if (text.labels) geom_text_repel else geom_label_repel + } else { + if (text.labels) geom_text else geom_label + } + + # Common arguments shared by all text geometries + label_args <- list( + data = vector_data, + mapping = aes( + x = .data[[xvar]], + y = .data[[yvar]], + label = .data[["vector_label"]] + ), + colour = label.color, + size = label.size, + inherit.aes = FALSE + ) + + # Add geometry-specific arguments + if (repel.labels) { + label_args <- c(label_args, list( + min.segment.length = min.segment.length, + box.padding = box.padding, + point.padding = point.padding, + force = force, + force_pull = force_pull, + max.time = max.time, + max.iter = max.iter, + max.overlaps = max.overlaps, + direction = direction, + seed = seed + )) + } + + # Add feature labels + p <- p + do.call(FUN, label_args) + + return(p) +} diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 7d914b50..ee788b7f 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -17,20 +17,231 @@ #' \code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}} #' object. #' +#' @param dimred \code{character scalar}. Name of the reduced dimension result +#' stored in \code{reducedDim(x)} to visualize. +#' #' @param ... Additional parameters for plotting. #' \itemize{ +#' \item \code{ncomponents}: \code{integer vector} of length 2 or +#' \code{integer scalar}. Specifies which ordination components are plotted. +#' If a scalar is provided, the first two components are used. +#' (Default: \code{2L}) +#' #' \item \code{colour.by}: \code{NULL} or \code{character scalar}. Specifies a -#' variable from \code{colData(x)} or \code{rowData(x)} which is used to -#' colour observations. (Default: \code{NULL}) +#' variable from \code{colData(x)} or a feature from \code{rownames(x)} used +#' to colour observations. Feature abundances are taken from +#' \code{assay.type}. (Default: \code{NULL}) +#' +#' \item \code{fill.by}: \code{NULL} or \code{character scalar}. Specifies a +#' variable from \code{colData(x)} or a feature from \code{rownames(x)} used +#' to fill observations or ellipses. Feature abundances are taken from +#' \code{assay.type}. Cannot be used together with +#' \code{add.density = TRUE}. (Default: \code{NULL}) +#' +#' \item \code{shape.by}: \code{NULL} or \code{character scalar}. Specifies a +#' categorical variable from \code{colData(x)} used for point shapes. +#' (Default: \code{NULL}) +#' +#' \item \code{size.by}: \code{NULL} or \code{character scalar}. Specifies a +#' variable from \code{colData(x)} used for point sizes. +#' (Default: \code{NULL}) +#' +#' \item \code{group.by}: \code{NULL} or \code{character scalar}. Specifies a +#' categorical variable from \code{colData(x)} used for grouping when drawing +#' ellipses, centroids and centroid vectors. (Default: \code{NULL}) +#' +#' \item \code{linetype.by}: \code{NULL} or \code{character scalar}. +#' Specifies a categorical variable from \code{colData(x)} used for ellipse +#' line types. (Default: \code{NULL}) +#' +#' \item \code{pair.by}: \code{NULL} or \code{character scalar}. Specifies a +#' variable from \code{colData(x)} identifying observations that should be +#' connected by lines. (Default: \code{NULL}) +#' +#' \item \code{sort.by}: \code{NULL} or \code{character scalar}. Specifies a +#' variable from \code{colData(x)} used to order observations before drawing +#' connecting lines. (Default: \code{NULL}) +#' +#' \item \code{facet.by}: \code{NULL} or \code{character scalar}. Specifies a +#' categorical variable from \code{colData(x)} used to split the plot into +#' facets. (Default: \code{NULL}) +#' +#' \item \code{assay.type}: \code{character scalar}. Name of the assay used +#' when \code{colour.by} or \code{fill.by} specifies a feature. +#' (Default: \code{"counts"}) +#' +#' \item \code{add.points}: \code{logical scalar}. Whether to draw sample +#' points. (Default: \code{TRUE}) +#' +#' \item \code{add.ellipse}: \code{logical scalar}. Whether to draw confidence +#' ellipses around groups. (Default: \code{FALSE}) +#' +#' \item \code{add.density}: \code{logical scalar}. Whether to draw a +#' two-dimensional density estimate in the background. +#' (Default: \code{FALSE}) +#' +#' \item \code{add.centroids}: \code{logical scalar}. Whether to draw group +#' centroids. (Default: \code{FALSE}) +#' +#' \item \code{add.centroids.lines}: \code{logical scalar}. Whether to connect +#' observations to their group centroids. (Default: \code{FALSE}) +#' +#' \item \code{add.vectors}: \code{logical scalar}. Whether to draw vectors +#' from the global centroid to group centroids. +#' (Default: \code{FALSE}) +#' +#' \item \code{add.rotation}: \code{logical scalar}. Whether to draw rotation +#' (species score) coordinates if available in the ordination result. +#' (Default: \code{add.species}) +#' +#' \item \code{add.species}: \code{logical scalar}. Alias for +#' \code{add.rotation}. (Default: \code{FALSE}) +#' +#' \item \code{add.expl.var}: \code{logical scalar}. Whether to append the +#' percentage of explained variance to the axis labels when available. +#' (Default: \code{FALSE}) +#' +#' \item \code{scales}: \code{character scalar}. Scaling used for faceted +#' plots. Passed to \code{ggplot2::facet_wrap()}. (Default: +#' \code{"fixed"}) +#' +#' \item \code{xlab}, \code{ylab}: \code{character scalar}. Axis labels. +#' Defaults to the ordination component names. +#' +#' \item \code{panel.by.eigen}: \code{logical scalar}. Whether to scale the +#' panel aspect ratio according to the eigenvalues of the ordination when +#' available. (Default: \code{TRUE}) +#' +#' \item \code{point.shape}: Shape used for points. +#' (Default: \code{19}) +#' +#' \item \code{point.alpha}: \code{numeric scalar}. Transparency of points. +#' Must be between 0 and 1. (Default: \code{0.4}) +#' +#' \item \code{ellipse.alpha}: \code{numeric scalar}. Transparency of ellipse +#' fills. Must be between 0 and 1. (Default: \code{0.2}) +#' +#' \item \code{ellipse.linewidth}: \code{numeric scalar}. Line width of +#' ellipse borders. (Default: \code{0.5} or \code{0} when +#' \code{fill.by} is specified.) +#' +#' \item \code{ellipse.linetype}: Integer specifying the ellipse line type. +#' (Default: \code{1}) +#' +#' \item \code{confidence.level}: \code{numeric scalar}. Confidence level used +#' for ellipse calculation. Must be between 0 and 1. +#' (Default: \code{0.95}) +#' +#' \item \code{adjust}: \code{numeric scalar}. Multiplicative adjustment for +#' the bandwidth used in the background density estimate. +#' (Default: \code{1}) #' } #' #' @examples #' data("Tito2024QMP") #' tse <- Tito2024QMP #' +#' # Compute relative abundances and an MDS ordination +#' tse <- transformAssay(tse, method = "relabundance") +#' tse <- addMDS(tse, assay.type = "relabundance", method = "bray", ncomponents = 50) +#' +#' # Basic ordination plot +#' plotOrdination(tse, "MDS") +#' +#' # Colour samples by a sample-level variable +#' plotOrdination(tse, "MDS", colour.by = "diagnosis") +#' +#' # Colour samples by the abundance of a single feature +#' plotOrdination( +#' tse, "MDS", +#' colour.by = rownames(tse)[1], +#' assay.type = "relabundance") +#' +#' # Use multiple aesthetics simultaneously +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' shape.by = "colonoscopy" +#' ) +#' +#' # Add confidence ellipses +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' group.by = "diagnosis", +#' add.ellipse = TRUE +#' ) +#' +#' # Fill ellipses instead of colouring them +#' plotOrdination( +#' tse, "MDS", +#' fill.by = "diagnosis", +#' group.by = "diagnosis", +#' add.ellipse = TRUE +#' ) +#' +#' # Show explained variance in the axis labels +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' add.expl.var = TRUE +#' ) +#' +#' # Add group centroids +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' group.by = "diagnosis", +#' add.centroids = TRUE +#' ) +#' +#' # Connect samples to their group centroids +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' group.by = "diagnosis", +#' add.centroids.lines = TRUE +#' ) +#' +#' # Show vectors from the global centroid to each group centroid +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' group.by = "diagnosis", +#' add.vectors = TRUE +#' ) +#' +#' # Draw a background density estimate +#' plotOrdination( +#' tse, "MDS", +#' add.density = TRUE +#' ) +#' +#' # Split the plot into facets +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' facet.by = "colonoscopy" +#' ) +#' +#' # Plot different ordination components +#' plotOrdination( +#' tse, "MDS", +#' ncomponents = c(2, 3) +#' ) +#' +#' # Customize point appearance +#' plotOrdination( +#' tse, "MDS", +#' colour.by = "diagnosis", +#' point.shape = 17, +#' point.alpha = 0.8 +#' ) +#' #' @seealso #' \itemize{ #' \item \code{\link[scater:plotReducedDim]{scater::plotReducedDim}} +#' \item \code{\link[=plotCCA]{plotCCA}} #' } #' NULL @@ -64,21 +275,10 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), assay.type = "counts", add.points = TRUE, add.ellipse = FALSE, add.density = FALSE, add.centroids = FALSE, add.centroids.lines = FALSE, add.vectors = FALSE, - add.rotation = FALSE, add.expl.var = FALSE, + add.rotation = add.species, add.species = FALSE, add.expl.var = FALSE, ...){ - # Check if there are any reduced dim present - if( length(reducedDims(x)) == 0L ){ - stop("No data present in reducedDim(x).", call. = FALSE) - } # Check that dimred can be found - is_name <- .is_a_string(dimred) && dimred %in% reducedDimNames(x) - is_index <- .is_an_integer(dimred) && dimred > 0L && - dimred <= length(reducedDims(x)) - if( !(is_name || is_index) ){ - stop("'dimred' must specify data from reducedDim(x). It must be one ", - "of the following options: '", - paste0(reducedDimNames(x), collapse = "', '"), "'", call. = FALSE) - } + .check_dimred_present(dimred, x) # Check that ncomponents is correct. We can only visualize 2 components. if( .is_an_integer(ncomponents) ){ ncomponents <- seq_len(ncomponents) @@ -91,16 +291,15 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Check aesthetic variables - temp <- .check_metadata_variable(tse, colour.by, FALSE, TRUE, FALSE, TRUE) - temp <- .check_metadata_variable(tse, fill.by, FALSE, TRUE, FALSE, TRUE) - temp <- .check_metadata_variable(tse, shape.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, size.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable( - tse, linetype.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, group.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, pair.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, sort.by, FALSE, TRUE, FALSE, FALSE) - temp <- .check_metadata_variable(tse, facet.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, colour.by, FALSE, TRUE, FALSE, TRUE) + temp <- .check_metadata_variable(x, fill.by, FALSE, TRUE, FALSE, TRUE) + temp <- .check_metadata_variable(x, shape.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, size.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, linetype.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, group.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, pair.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, sort.by, FALSE, TRUE, FALSE, FALSE) + temp <- .check_metadata_variable(x, facet.by, FALSE, TRUE, FALSE, FALSE) # If colour.by specifies rowname, we check assay.type as the abundance # values are used for coloring if( !is.null(colour.by) && colour.by %in% rownames(x) ){ @@ -165,7 +364,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), add.expl.var = FALSE, ...){ # Get data and store the original attributes that might include rotation # data, for instance - df <- reducedDim(x, dimred) + df <- reducedDim(x, dimred)[, ncomponents] orig_attributes <- attributes(df) orig_attributes <- orig_attributes[ !names(orig_attributes) %in% c("dim", "dimnames") ] @@ -176,25 +375,24 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } df <- df |> as.data.frame() # Take only 2 specified columns - df <- df[, ncomponents] x_var <- colnames(df)[[1L]] y_var <- colnames(df)[[2L]] - # Get rotation data adnd put it in correct format + # Get rotation data and put it in correct format rotation <- NULL - rotation_names <- c("rotation") + rotation_names <- c("rotation", "species") if( any(rotation_names %in% names(orig_attributes)) ){ - rotation_names <- rotation_names[[1L]] + rotation_names <- rotation_names[ + rotation_names %in% names(orig_attributes)][[1L]] rotation <- orig_attributes[[rotation_names]] |> as.data.frame() rotation <- rotation[, ncomponents, drop = FALSE] colnames(rotation) <- colnames(df) } - expl_var_name <- c("eig") xlab <- x_var ylab <- y_var - if( add.expl.var && any(expl_var_name %in% names(orig_attributes)) ){ - eigen <- orig_attributes[expl_var_name][[1L]] + eigen <- orig_attributes[["eig"]] + if( add.expl.var && !is.null(eigen) ){ xlab <- paste0( xlab, " (", round(eigen[ncomponents][[1L]], 1), @@ -284,6 +482,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), ) attr(df, "rotation") <- rotation attr(df, "centroids") <- df_centroids + attr(df, "eigen") <- eigen return(df) } @@ -291,7 +490,8 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), .ordination_plotter <- function( df, scales = "fixed", add.points = TRUE, add.ellipse = FALSE, add.density = FALSE, add.centroids = FALSE, add.centroids.lines = FALSE, - add.vectors = FALSE, add.rotation = FALSE, ...){ + add.vectors = FALSE, add.rotation = add.species, add.species = FALSE, + ...){ # Initialize the plot p <- ggplot(df, aes( x = .data[[attributes(df)[["x"]]]], @@ -345,11 +545,6 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Adjust theme p <- .adjust_ordination_theme(p, df, ...) - - # Enforce same scale to x and y axis. Without equal scale, the results and - # interpretations might be misleading - p <- p + coord_equal() - return(p) } @@ -473,6 +668,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # This methods creates vectors that start from global mean and ends to group # centroids. This shows how the covariate correlates with the ordination. +#' @importFrom ggrepel geom_label_repel .add_centroids_vector <- function(p, df, grouping_var, ...){ df_centroids <- attributes(df)[["centroids"]] # Visualize vectors @@ -486,7 +682,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), size = 1 ) # Add label to denote which vector belongs to which group - p <- p + ggrepel::geom_label_repel( + p <- p + geom_label_repel( data = df_centroids, mapping = aes( x = x_centroid, @@ -500,7 +696,6 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # ordination. .add_rotation <- function(p, df){ # The points are added only if the rotation is present in the data - rot_names <- c("rotation") if( !is.null(attributes(df)[["rotation"]]) ){ df_species <- attributes(df)[["rotation"]] # Add points @@ -549,6 +744,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), p, df, xlab = attributes(df)[["xlab"]], ylab = attributes(df)[["ylab"]], + panel.by.eigen = TRUE, ...){ if( !.is_a_string(xlab) ){ stop("'xlab' must be a single character value.", call. = FALSE) @@ -556,6 +752,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !.is_a_string(ylab) ){ stop("'ylab' must be a single character value.", call. = FALSE) } + if( !.is_a_bool(panel.by.eigen) ){ + stop("'panel.by.eigen' must be TRUE or FALSE.", call. = FALSE) + } # p <- p + theme_classic() p <- p + labs(x = xlab, y = ylab) @@ -574,5 +773,28 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !is.null(attributes(df)[["linetype.by"]]) ){ p <- p + labs(linetype = attributes(df)[["linetype.by"]]) } + + # Use equal scaling on both axes so that one unit on the x-axis has the same + # physical length as one unit on the y-axis. This preserves distances and + # angles in the ordination and avoids visual distortion. + p <- p + coord_equal() + # Make the panel dimensions proportional to the eigenvalues so that the + # relative lengths of the axes reflect the variation explained by each + # ordination axis. + if( panel.by.eigen && !is.null(attributes(df)[["eigen"]]) ){ + eig <- attributes(df)[["eigen"]] + p <- p + theme(aspect.ratio = eig[2] / eig[1]) + } + + # Adjust colors + name <- if(!is.null(attributes(df)[["colour.by"]])) + attributes(df)[["colour.by"]] else attributes(df)[["fill.by"]] + vals <- if(!is.null(attributes(df)[["colour.by"]])) + df[[attributes(df)[["colour.by"]]]] else attributes(df)[["fill.by"]] + p <- .resolve_plot_colours( + p, vals, name, + fill = !is.null(attributes(df)[["fill.by"]]) + ) + return(p) } diff --git a/R/utils.R b/R/utils.R index 75ee2456..1e41bf41 100644 --- a/R/utils.R +++ b/R/utils.R @@ -24,6 +24,8 @@ TAXONOMY_RANKS <- mia:::TAXONOMY_RANKS .check_rowTree_present <- mia:::.check_rowTree_present .check_colTree_present <- mia:::.check_colTree_present .merge_features <- mia:::.merge_features +.check_dimred_present <- mia:::.check_dimred_present +.check_metadata_present <- mia:::.check_metadata_present .norm_label <- function(label, x){ if(!is.null(label)){ diff --git a/man/plotCCA.Rd b/man/plotCCA.Rd index 03377eb0..e8021ddf 100644 --- a/man/plotCCA.Rd +++ b/man/plotCCA.Rd @@ -4,9 +4,7 @@ \alias{plotCCA} \alias{plotRDA} \alias{plotCCA,SingleCellExperiment-method} -\alias{plotCCA,matrix-method} \alias{plotRDA,SingleCellExperiment-method} -\alias{plotRDA,matrix-method} \title{Plot RDA or CCA object} \usage{ plotCCA(x, ...) @@ -15,17 +13,12 @@ plotRDA(x, ...) \S4method{plotCCA}{SingleCellExperiment}(x, dimred, ...) -\S4method{plotCCA}{matrix}(x, ...) - \S4method{plotRDA}{SingleCellExperiment}(x, dimred, ...) - -\S4method{plotRDA}{matrix}(x, ...) } \arguments{ \item{x}{a \code{\link[TreeSummarizedExperiment:TreeSummarizedExperiment-constructor]{TreeSummarizedExperiment}} -or a matrix of weights. The latter is returned as output from -\code{\link[mia:runCCA]{getRDA}}.} +object.} \item{...}{additional parameters for plotting, inherited from \code{\link[scater:plotReducedDim]{plotReducedDim}}, @@ -140,14 +133,17 @@ tse <- addRDA( formula = assay ~ ClinicalStatus + Gender + Age, distance = "bray", na.action = na.exclude - ) +) suppressWarnings({ # Create RDA plot coloured by variable plotRDA(tse, "RDA", colour.by = "ClinicalStatus") -# Create RDA plot with empty ellipses -plotRDA(tse, "RDA", colour.by = "ClinicalStatus", add.ellipse = "colour") +# Create RDA plot with ellipses +plotRDA( + tse, "RDA", colour.by = "ClinicalStatus", fill.by = "ClinicalStatus", + add.ellipse = TRUE +) # Create RDA plot with text encased in labels plotRDA(tse, "RDA", colour.by = "ClinicalStatus", vec.text = FALSE) @@ -157,18 +153,11 @@ plotRDA(tse, "RDA", colour.by = "ClinicalStatus", repel.labels = FALSE) # Create RDA plot without vectors plotRDA(tse, "RDA", colour.by = "ClinicalStatus", add.vectors = FALSE) - -# Calculate RDA as a separate object -rda_mat <- getRDA( - tse, - assay.type = "relabundance", - formula = assay ~ ClinicalStatus + Gender + Age, - distance = "bray", - na.action = na.exclude - ) - -# Create RDA plot from RDA matrix -plotRDA(rda_mat) }) } +\seealso{ +\itemize{ +\item \code{\link[=plotOrdination]{plotOrdination}} +} +} diff --git a/man/plotJointRPCA.Rd b/man/plotJointRPCA.Rd new file mode 100644 index 00000000..f0adbc7a --- /dev/null +++ b/man/plotJointRPCA.Rd @@ -0,0 +1,121 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/AllGenerics.R, R/plotJointRPCA.R +\name{plotJointRPCA} +\alias{plotJointRPCA} +\alias{plotJointRPCA,MultiAssayExperiment-method} +\alias{plotJointRPCA,SingleCellExperiment-method} +\title{Visualize Joint-RPCA results} +\usage{ +plotJointRPCA(x, ...) + +\S4method{plotJointRPCA}{MultiAssayExperiment}(x, dimred, ...) + +\S4method{plotJointRPCA}{SingleCellExperiment}(x, dimred, ...) +} +\arguments{ +\item{x}{A +\code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}} +or +\code{\link[MultiAssayExperiment:MultiAssayExperiment-class]{MultiAssayExperiment}} +object containing Joint-RPCA results.} + +\item{...}{Additional arguments passed to +\code{\link[=plotOrdination]{plotOrdination}} and to the feature-vector +plotting. +Commonly used Joint-RPCA-specific arguments include: +\itemize{ +\item \code{add.vectors}: \code{logical(1)} or \code{character}. If +\code{TRUE}, feature loading vectors are added. Alternatively, a character +vector can be supplied to display only features whose names match the given +pattern(s). (Default: \code{TRUE}) +\item \code{ntop}: \code{integer(1)} or \code{NULL}. Maximum number of +loading vectors shown per data layer. The vectors with the largest lengths +are retained. Set to \code{NULL} to display all vectors. (Default: \code{10}) +}} + +\item{dimred}{\code{character(1)} specifying the name of the Joint-RPCA result to plot.} +} +\value{ +A \code{ggplot2} object. +} +\description{ +Creates a two-dimensional ordination plot from Joint-RPCA results stored in a +\code{SingleCellExperiment} or \code{MultiAssayExperiment}. In addition to the +sample ordination, feature loadings can be visualized as vectors. +} +\details{ +This function is a wrapper around \code{\link[=plotOrdination]{plotOrdination}} +for Joint-RPCA results. Consequently, most graphical parameters are passed +directly to \code{plotOrdination()}, including options for colouring, +grouping, faceting and adding ellipses. See +\code{\link[=plotOrdination]{plotOrdination}} for a complete description of +these arguments. + +Feature loadings are plotted as vectors originating from the origin. By +default, only the longest loading vectors are shown for each data layer to +improve readability. +} +\examples{ +data("HintikkaXOData") +mae <- HintikkaXOData + +mae[[1]] <- transformAssay( + mae[[1]], + assay.type = "counts", + method = "rclr", + impute = FALSE +) +mae[[2]] <- transformAssay( + mae[[2]], + assay.type = "nmr", + method = "log10" +) + +mae <- addJointRPCA( + mae, + experiments = c(1, 2), + assay.types = c("rclr", "log10") +) + +# Basic Joint-RPCA plot +plotJointRPCA(mae, "JointRPCA") + +# Colour samples by metadata +plotJointRPCA(mae, "JointRPCA", colour.by = "Fat") + +# Add confidence ellipses +plotJointRPCA( + mae, + "JointRPCA", + colour.by = "Fat", + add.ellipse = TRUE +) + +# Show all loading vectors +plotJointRPCA( + mae, + "JointRPCA", + ntop = NULL +) + +## Display only selected feature vectors +plotJointRPCA( + mae, + "JointRPCA", + add.vectors = c("Bacteroides", "Roseburia") +) + +## Use boxed labels without repelling +plotJointRPCA( + mae, + "JointRPCA", + text.labels = FALSE, + repel.labels = FALSE +) + +} +\seealso{ +\itemize{ +\item \code{\link[=plotOrdination]{plotOrdination}} +} +} diff --git a/man/plotOrdination.Rd b/man/plotOrdination.Rd new file mode 100644 index 00000000..7be443e9 --- /dev/null +++ b/man/plotOrdination.Rd @@ -0,0 +1,253 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/AllGenerics.R, R/plotOrdination.R +\name{plotOrdination} +\alias{plotOrdination} +\alias{plotOrdination,SingleCellExperiment-method} +\title{Create ordination plot} +\usage{ +plotOrdination(x, ...) + +\S4method{plotOrdination}{SingleCellExperiment}(x, dimred, colour.by = color.by, color.by = NULL, ...) +} +\arguments{ +\item{x}{a +\code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}} +object.} + +\item{...}{Additional parameters for plotting. +\itemize{ +\item \code{ncomponents}: \code{integer vector} of length 2 or +\code{integer scalar}. Specifies which ordination components are plotted. +If a scalar is provided, the first two components are used. +(Default: \code{2L}) + +\item \code{colour.by}: \code{NULL} or \code{character scalar}. Specifies a +variable from \code{colData(x)} or a feature from \code{rownames(x)} used +to colour observations. Feature abundances are taken from +\code{assay.type}. (Default: \code{NULL}) + +\item \code{fill.by}: \code{NULL} or \code{character scalar}. Specifies a +variable from \code{colData(x)} or a feature from \code{rownames(x)} used +to fill observations or ellipses. Feature abundances are taken from +\code{assay.type}. Cannot be used together with +\code{add.density = TRUE}. (Default: \code{NULL}) + +\item \code{shape.by}: \code{NULL} or \code{character scalar}. Specifies a +categorical variable from \code{colData(x)} used for point shapes. +(Default: \code{NULL}) + +\item \code{size.by}: \code{NULL} or \code{character scalar}. Specifies a +variable from \code{colData(x)} used for point sizes. +(Default: \code{NULL}) + +\item \code{group.by}: \code{NULL} or \code{character scalar}. Specifies a +categorical variable from \code{colData(x)} used for grouping when drawing +ellipses, centroids and centroid vectors. (Default: \code{NULL}) + +\item \code{linetype.by}: \code{NULL} or \code{character scalar}. +Specifies a categorical variable from \code{colData(x)} used for ellipse +line types. (Default: \code{NULL}) + +\item \code{pair.by}: \code{NULL} or \code{character scalar}. Specifies a +variable from \code{colData(x)} identifying observations that should be +connected by lines. (Default: \code{NULL}) + +\item \code{sort.by}: \code{NULL} or \code{character scalar}. Specifies a +variable from \code{colData(x)} used to order observations before drawing +connecting lines. (Default: \code{NULL}) + +\item \code{facet.by}: \code{NULL} or \code{character scalar}. Specifies a +categorical variable from \code{colData(x)} used to split the plot into +facets. (Default: \code{NULL}) + +\item \code{assay.type}: \code{character scalar}. Name of the assay used +when \code{colour.by} or \code{fill.by} specifies a feature. +(Default: \code{"counts"}) + +\item \code{add.points}: \code{logical scalar}. Whether to draw sample +points. (Default: \code{TRUE}) + +\item \code{add.ellipse}: \code{logical scalar}. Whether to draw confidence +ellipses around groups. (Default: \code{FALSE}) + +\item \code{add.density}: \code{logical scalar}. Whether to draw a +two-dimensional density estimate in the background. +(Default: \code{FALSE}) + +\item \code{add.centroids}: \code{logical scalar}. Whether to draw group +centroids. (Default: \code{FALSE}) + +\item \code{add.centroids.lines}: \code{logical scalar}. Whether to connect +observations to their group centroids. (Default: \code{FALSE}) + +\item \code{add.vectors}: \code{logical scalar}. Whether to draw vectors +from the global centroid to group centroids. +(Default: \code{FALSE}) + +\item \code{add.rotation}: \code{logical scalar}. Whether to draw rotation +(species score) coordinates if available in the ordination result. +(Default: \code{add.species}) + +\item \code{add.species}: \code{logical scalar}. Alias for +\code{add.rotation}. (Default: \code{FALSE}) + +\item \code{add.expl.var}: \code{logical scalar}. Whether to append the +percentage of explained variance to the axis labels when available. +(Default: \code{FALSE}) + +\item \code{scales}: \code{character scalar}. Scaling used for faceted +plots. Passed to \code{ggplot2::facet_wrap()}. (Default: +\code{"fixed"}) + +\item \code{xlab}, \code{ylab}: \code{character scalar}. Axis labels. +Defaults to the ordination component names. + +\item \code{panel.by.eigen}: \code{logical scalar}. Whether to scale the +panel aspect ratio according to the eigenvalues of the ordination when +available. (Default: \code{TRUE}) + +\item \code{point.shape}: Shape used for points. +(Default: \code{19}) + +\item \code{point.alpha}: \code{numeric scalar}. Transparency of points. +Must be between 0 and 1. (Default: \code{0.4}) + +\item \code{ellipse.alpha}: \code{numeric scalar}. Transparency of ellipse +fills. Must be between 0 and 1. (Default: \code{0.2}) + +\item \code{ellipse.linewidth}: \code{numeric scalar}. Line width of +ellipse borders. (Default: \code{0.5} or \code{0} when +\code{fill.by} is specified.) + +\item \code{ellipse.linetype}: Integer specifying the ellipse line type. +(Default: \code{1}) + +\item \code{confidence.level}: \code{numeric scalar}. Confidence level used +for ellipse calculation. Must be between 0 and 1. +(Default: \code{0.95}) + +\item \code{adjust}: \code{numeric scalar}. Multiplicative adjustment for +the bandwidth used in the background density estimate. +(Default: \code{1}) +}} + +\item{dimred}{\code{character scalar}. Name of the reduced dimension result +stored in \code{reducedDim(x)} to visualize.} +} +\value{ +A \code{ggplot2} object. +} +\description{ +Ordinaton plotter +} +\details{ +Creates ordination plot +} +\examples{ +data("Tito2024QMP") +tse <- Tito2024QMP + +# Compute relative abundances and an MDS ordination +tse <- transformAssay(tse, method = "relabundance") +tse <- addMDS(tse, assay.type = "relabundance", method = "bray", ncomponents = 50) + +# Basic ordination plot +plotOrdination(tse, "MDS") + +# Colour samples by a sample-level variable +plotOrdination(tse, "MDS", colour.by = "diagnosis") + +# Colour samples by the abundance of a single feature +plotOrdination( + tse, "MDS", + colour.by = rownames(tse)[1], + assay.type = "relabundance") + +# Use multiple aesthetics simultaneously +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + shape.by = "colonoscopy" +) + +# Add confidence ellipses +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + group.by = "diagnosis", + add.ellipse = TRUE +) + +# Fill ellipses instead of colouring them +plotOrdination( + tse, "MDS", + fill.by = "diagnosis", + group.by = "diagnosis", + add.ellipse = TRUE +) + +# Show explained variance in the axis labels +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + add.expl.var = TRUE +) + +# Add group centroids +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + group.by = "diagnosis", + add.centroids = TRUE +) + +# Connect samples to their group centroids +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + group.by = "diagnosis", + add.centroids.lines = TRUE +) + +# Show vectors from the global centroid to each group centroid +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + group.by = "diagnosis", + add.vectors = TRUE +) + +# Draw a background density estimate +plotOrdination( + tse, "MDS", + add.density = TRUE +) + +# Split the plot into facets +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + facet.by = "colonoscopy" +) + +# Plot different ordination components +plotOrdination( + tse, "MDS", + ncomponents = c(2, 3) +) + +# Customize point appearance +plotOrdination( + tse, "MDS", + colour.by = "diagnosis", + point.shape = 17, + point.alpha = 0.8 +) + +} +\seealso{ +\itemize{ +\item \code{\link[scater:plotReducedDim]{scater::plotReducedDim}} +\item \code{\link[=plotCCA]{plotCCA}} +} +} diff --git a/tests/testthat/test-plotJointRPCA.R b/tests/testthat/test-plotJointRPCA.R new file mode 100644 index 00000000..88100604 --- /dev/null +++ b/tests/testthat/test-plotJointRPCA.R @@ -0,0 +1,287 @@ + +.make_joint <- function() { + + tse <- makeTSE(nrow = 8, ncol = 4) + assayNames(tse) <- "counts" + + rd <- matrix(rnorm(8), ncol = 2) + class(rd) <- c("JointRPCA", class(rd)) + + rotation <- matrix(rnorm(16), ncol = 2) + rownames(rotation) <- paste0("feature", seq_len(nrow(rotation))) + + attr(rd, "rotation") <- rotation + attr(rd, "n_features") <- c(Layer1 = 4, Layer2 = 4) + + reducedDim(tse, "JointRPCA") <- rd + + return(tse) +} + +test_that("returns ggplot", { + + p <- plotJointRPCA(.make_joint(), "JointRPCA") + + expect_s3_class(p, "ggplot") + +}) + +test_that("FALSE disables vectors", { + + expect_no_error( + + plotJointRPCA( + .make_joint(), + "JointRPCA", + add.vectors = FALSE + ) + + ) + +}) + +test_that("character vectors accepted", { + + expect_no_error( + + plotJointRPCA( + .make_joint(), + "JointRPCA", + add.vectors = "feature" + ) + + ) + +}) + +test_that("invalid add.vectors errors", { + + expect_error( + + plotJointRPCA( + .make_joint(), + "JointRPCA", + add.vectors = 1 + ), + + "add.vectors" + + ) + +}) + +test_that("non JointRPCA reducedDim errors", { + + tse <- makeTSE(nrow = 5, ncol = 2) + + reducedDim(tse,"PCA") <- matrix(rnorm(4),2) + + expect_error( + + plotJointRPCA( + tse, + "PCA" + ), + + "JointRPCA" + + ) + +}) + +test_that("returns data.frame", { + + x <- .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + TRUE + ) + + expect_s3_class(x,"data.frame") + +}) + +test_that("FALSE returns NULL", { + + expect_null( + + .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + FALSE + ) + + ) + +}) + +test_that("feature filtering works", { + + x <- .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + "feature1" + ) + + expect_true( + + all(grepl( + "feature1", + x$vector_label + )) + + ) + +}) + +test_that("case insensitive filtering works", { + + x <- .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + "FEATURE", + ignore.case = TRUE + ) + + expect_gt(nrow(x),0) + +}) + +test_that("ntop limits each layer", { + + x <- .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + TRUE, + ntop = 2 + ) + + expect_true( + + all(table(x$Layer) <= 2) + + ) + +}) + +test_that("ntop NULL keeps all", { + + x <- .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + TRUE, + ntop = NULL + ) + + expect_equal( + + nrow(x), + + 8 + + ) + +}) + +test_that("ignore.case validated", { + + expect_error( + + .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + TRUE, + ignore.case = 1 + ), + + "ignore.case" + + ) + +}) + +test_that("ntop validated", { + + expect_error( + + .get_joint_rpca_vector_data( + .make_joint(), + "JointRPCA", + TRUE, + ntop = "a" + ), + + "ntop" + + ) + +}) + +.make_vectors <- function() { + + data.frame( + PC1 = c(.2,.5), + PC2 = c(.4,.1), + Layer = c("A","B"), + vector_label = c("x","y") + ) + +} + +test_that("returns ggplot", { + + p <- ggplot() + + p2 <- .add_ordination_vectors( + p, + .make_vectors() + ) + + expect_s3_class(p2,"ggplot") + + expect_length( + + p2$layers, + + 2 + + ) + +}) + +test_that("supports geom_label", { + + p <- .add_ordination_vectors( + ggplot(), + .make_vectors(), + text.labels = FALSE + ) + + expect_s3_class(p,"ggplot") + +}) + +test_that("supports non-repel labels", { + + p <- .add_ordination_vectors( + ggplot(), + .make_vectors(), + repel.labels = FALSE + ) + + expect_s3_class(p,"ggplot") + +}) + +test_that("supports label without repel", { + + p <- .add_ordination_vectors( + ggplot(), + .make_vectors(), + repel.labels = FALSE, + text.labels = FALSE + ) + + expect_s3_class(p,"ggplot") + +}) diff --git a/tests/testthat/test-plotOrdination.R b/tests/testthat/test-plotOrdination.R new file mode 100644 index 00000000..4b85796a --- /dev/null +++ b/tests/testthat/test-plotOrdination.R @@ -0,0 +1,228 @@ + +.make_tse <- function() { + tse <- makeTSE(ncol = 2, nrow = 10) + assayNames(tse) <- "counts" + + reducedDim(tse, "PCA") <- matrix( + seq_len(6), + ncol = 3, + dimnames = list(colnames(tse), c("PC1", "PC2", "PC3")) + ) + + attr(reducedDim(tse, "PCA"), "eig") <- c(60, 30, 10) + + colData(tse)$numeric <- c(1, 2) + colData(tse)$time <- c(1, 2) + + return(tse) +} + +test_that("defaults are returned", { + + tse <- .make_tse() + + args <- .check_ordination_input(tse, "PCA") + + expect_identical(args$x, tse) + expect_identical(args$dimred, "PCA") + expect_identical(args$ncomponents, 1:2) + expect_true(args$add.expl.var == FALSE) + +}) + +test_that("scalar ncomponents expands", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = 2L + ) + + expect_identical(args$ncomponents, c(1L, 2L)) + +}) + +test_that("vector ncomponents is preserved", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = c(2L, 3L) + ) + + expect_identical(args$ncomponents, c(2L, 3L)) + +}) + +test_that("three components are rejected", { + + expect_error( + .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = 1:3 + ), + "ncomponents" + ) + +}) + +test_that("component outside range is rejected", { + + expect_error( + .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = c(1L, 4L) + ), + "ncomponents" + ) + +}) + +test_that("non-integer components are rejected", { + + expect_error( + .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = c(1, 2.5) + ), + "ncomponents" + ) + +}) + +test_that("density and fill are mutually exclusive", { + + expect_error( + .check_ordination_input( + .make_tse(), + "PCA", + fill.by = "group", + add.density = TRUE + ), + "Both 'add.density' and 'fill.by' cannot be specified simultaneously." + ) + +}) + +for (arg in c( + "add.points", + "add.ellipse", + "add.density", + "add.centroids", + "add.centroids.lines", + "add.vectors", + "add.rotation", + "add.expl.var" +)) { + + test_that(paste(arg, "must be logical"), { + + x <- list( + x = .make_tse(), + dimred = "PCA" + ) + + x[[arg]] <- 1 + + expect_error( + do.call(.check_ordination_input, x), + arg + ) + + }) + +} + +test_that("extra arguments are propagated", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + point.alpha = .5, + point.shape = 17 + ) + + expect_identical(args$point.alpha, .5) + expect_identical(args$point.shape, 17) + +}) + +test_that("metadata arguments are propagated", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + colour.by = "group", + fill.by = "group", + shape.by = "group", + size.by = "numeric", + group.by = "group", + linetype.by = "group", + pair.by = "ID", + sort.by = "time", + facet.by = "group" + ) + + expect_identical(args$colour.by, "group") + expect_identical(args$fill.by, "group") + expect_identical(args$shape.by, "group") + expect_identical(args$size.by, "numeric") + expect_identical(args$group.by, "group") + expect_identical(args$linetype.by, "group") + expect_identical(args$pair.by, "ID") + expect_identical(args$sort.by, "time") + expect_identical(args$facet.by, "group") + +}) + +test_that("plotOrdination returns a ggplot", { + p <- plotOrdination(.make_tse(), "PCA") + expect_s3_class(p, "ggplot") +}) + +test_that("integer ncomponents of length one expands to first two components", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = 2L + ) + + expect_identical(args$ncomponents, c(1L, 2L)) + +}) + +test_that("component order is preserved", { + + args <- .check_ordination_input( + .make_tse(), + "PCA", + ncomponents = c(3L, 1L) + ) + + expect_identical(args$ncomponents, c(3L, 1L)) + +}) + +test_that("all plotting flags can be enabled", { + + expect_no_error( + + .check_ordination_input( + .make_tse(), + "PCA", + add.points = TRUE, + add.ellipse = TRUE, + add.centroids = TRUE, + add.centroids.lines = TRUE, + add.vectors = TRUE, + add.rotation = TRUE, + add.expl.var = TRUE + ) + + ) + +}) From c1dfb7214a7833377d2b2838e8c654f98991c638 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 22:50:29 +0300 Subject: [PATCH 10/28] up --- DESCRIPTION | 2 +- R/plotCCA.R | 29 ----------------------------- 2 files changed, 1 insertion(+), 30 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 98c44e1e..630b0a7d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: miaViz Title: Microbiome Analysis Plotting and Visualization -Version: 1.21.1 +Version: 1.21.2 Authors@R: c(person(given = "Tuomas", family = "Borman", role = c("aut", "cre"), email = "tuomas.v.borman@utu.fi", diff --git a/R/plotCCA.R b/R/plotCCA.R index 1dc9fc3e..fdceaeec 100644 --- a/R/plotCCA.R +++ b/R/plotCCA.R @@ -206,35 +206,6 @@ setMethod("plotRDA", signature = c(x = "SingleCellExperiment"), return(object) } -# The data can include constrained and unconstrained axes. This function subsets -# the data so that it includes only constrained axes. -.subset_constrained_rda <- function(reduced_dim){ - # Get only the indices of constrained ones, i.e., first set of axes. - # The colnames are in format, constrained_axis1, ca2, ca3..., unconstrained - # axis1, uca2, ... - comp_num <- as.numeric(gsub("\\D", "", colnames(reduced_dim))) - ind <- which( cumsum(comp_num == 1) <= 1 ) - # If there were problems, it might be that the names are just arbitrary. - # Then take all the columns. - if( !(length(ind) > 0L && all(diff(ind) == 1L)) ){ - ind <- seq_len(ncol(reduced_dim)) - } - # Preserve attributes - attributes <- attributes(reduced_dim) - attributes <- attributes[ !names(attributes) %in% c("dim", "dimnames") ] - # Subset the data so that it includes only constrained axes - reduced_dim <- reduced_dim[ , ind, drop = FALSE] - if( "biplot" %in% names(attributes) ){ - attributes[["biplot"]] <- attributes[["biplot"]][ , ind, drop = FALSE] - } - if( "eig" %in% names(attributes) ){ - attributes[["eig"]] <- attributes[["eig"]][ind] - } - # Add attributes back - attributes(reduced_dim) <- c(attributes(reduced_dim), attributes) - return(reduced_dim) -} - # This function retrieves data for creating vectors. Moreover, it wrangles the # vector data and controls what information is added to vector text or labels. .get_rda_vector_data <- function( From 557ed6f7ffa1df49c4398ff68d5c0d2f837f65ae Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:05:36 +0300 Subject: [PATCH 11/28] up --- R/plotBoxplot.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/plotBoxplot.R b/R/plotBoxplot.R index 3b65cfb1..4e7d739d 100644 --- a/R/plotBoxplot.R +++ b/R/plotBoxplot.R @@ -1112,8 +1112,8 @@ setMethod("plotBoxplot", signature = c(object = "SummarizedExperiment"), x_point <- y_point <- NULL args <- list( mapping = aes( - x = .data[[x]], - y = .data[[y]], + x = x_point, + y = y_point, colour = if(!is.null(attributes(df)[["colour.by"]])) .data[[attributes(df)[["colour.by"]]]], shape = if(!is.null(attributes(df)[["shape.by"]])) From ef7cf5a912bbe6fbd44efc6aa6d8371d8be0d512 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:09:43 +0300 Subject: [PATCH 12/28] up --- R/plotBoxplot.R | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/R/plotBoxplot.R b/R/plotBoxplot.R index 4e7d739d..4630934a 100644 --- a/R/plotBoxplot.R +++ b/R/plotBoxplot.R @@ -1332,5 +1332,16 @@ setMethod("plotBoxplot", signature = c(object = "SummarizedExperiment"), if( !is.null(attributes(df)[["size.by"]]) ){ p <- p + labs(shape = attributes(df)[["size.by"]]) } + + # Adjust colors + name <- if(!is.null(attributes(df)[["colour.by"]])) + attributes(df)[["colour.by"]] else attributes(df)[["fill.by"]] + vals <- if(!is.null(attributes(df)[["colour.by"]])) + df[[attributes(df)[["colour.by"]]]] else attributes(df)[["fill.by"]] + p <- .resolve_plot_colours( + p, vals, name, + fill = !is.null(attributes(df)[["fill.by"]]) + ) + return(p) } From 8454640477b71e7bc866247dcf34aadf959a3e7d Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:23:17 +0300 Subject: [PATCH 13/28] up --- R/plotBoxplot.R | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/R/plotBoxplot.R b/R/plotBoxplot.R index 4630934a..95371c7e 100644 --- a/R/plotBoxplot.R +++ b/R/plotBoxplot.R @@ -1106,14 +1106,16 @@ setMethod("plotBoxplot", signature = c(object = "SummarizedExperiment"), # This function adds points to plot .add_points_layer <- function( - p, df, point.alpha = 0.65, point.size = 2, point.shape = 19L, + p, df, + x = NULL, y = NULL, + point.alpha = 0.65, point.size = 2, point.shape = 19L, point.colour = point.color, point.color = "grey70", ...){ # To disable "no visible binding for global variable" message in cmdcheck x_point <- y_point <- NULL args <- list( mapping = aes( - x = x_point, - y = y_point, + x = if(!is.null(x)) .data[[x]] else x_point, + y = if(!is.null(y)) .data[[y]] else y_point, colour = if(!is.null(attributes(df)[["colour.by"]])) .data[[attributes(df)[["colour.by"]]]], shape = if(!is.null(attributes(df)[["shape.by"]])) From a3c4a252cf8097ee27bfb59db79cb18152f26684 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:32:36 +0300 Subject: [PATCH 14/28] up --- R/plotBoxplot.R | 22 ++++++++++++++-------- R/plotOrdination.R | 22 ++++++++++++++-------- 2 files changed, 28 insertions(+), 16 deletions(-) diff --git a/R/plotBoxplot.R b/R/plotBoxplot.R index 95371c7e..89ab1739 100644 --- a/R/plotBoxplot.R +++ b/R/plotBoxplot.R @@ -1336,14 +1336,20 @@ setMethod("plotBoxplot", signature = c(object = "SummarizedExperiment"), } # Adjust colors - name <- if(!is.null(attributes(df)[["colour.by"]])) - attributes(df)[["colour.by"]] else attributes(df)[["fill.by"]] - vals <- if(!is.null(attributes(df)[["colour.by"]])) - df[[attributes(df)[["colour.by"]]]] else attributes(df)[["fill.by"]] - p <- .resolve_plot_colours( - p, vals, name, - fill = !is.null(attributes(df)[["fill.by"]]) - ) + if( !is.null(attributes(df)[["colour.by"]]) ){ + name <- attributes(df)[["colour.by"]] + vals <- df[[name]] + p <- .resolve_plot_colours( + p, vals, name, + ) + } + if( !is.null(attributes(df)[["fill.by"]]) ){ + name <- attributes(df)[["fill.by"]] + vals <- df[[name]] + p <- .resolve_plot_colours( + p, vals, name, fill = TRUE + ) + } return(p) } diff --git a/R/plotOrdination.R b/R/plotOrdination.R index ee788b7f..1bb21e38 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -787,14 +787,20 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), } # Adjust colors - name <- if(!is.null(attributes(df)[["colour.by"]])) - attributes(df)[["colour.by"]] else attributes(df)[["fill.by"]] - vals <- if(!is.null(attributes(df)[["colour.by"]])) - df[[attributes(df)[["colour.by"]]]] else attributes(df)[["fill.by"]] - p <- .resolve_plot_colours( - p, vals, name, - fill = !is.null(attributes(df)[["fill.by"]]) - ) + if( !is.null(attributes(df)[["colour.by"]]) ){ + name <- attributes(df)[["colour.by"]] + vals <- df[[name]] + p <- .resolve_plot_colours( + p, vals, name, + ) + } + if( !is.null(attributes(df)[["fill.by"]]) ){ + name <- attributes(df)[["fill.by"]] + vals <- df[[name]] + p <- .resolve_plot_colours( + p, vals, name, fill = TRUE + ) + } return(p) } From 82985452377c8a173f12ee7473651bbbd99fe559 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:57:02 +0300 Subject: [PATCH 15/28] up --- R/plotOrdination.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 1bb21e38..9f498d5f 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -364,11 +364,11 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), add.expl.var = FALSE, ...){ # Get data and store the original attributes that might include rotation # data, for instance - df <- reducedDim(x, dimred)[, ncomponents] + df <- reducedDim(x, dimred) orig_attributes <- attributes(df) orig_attributes <- orig_attributes[ !names(orig_attributes) %in% c("dim", "dimnames") ] - + df <- df[, ncomponents] # Add colnames if they are not present if( is.null(colnames(df)) ){ colnames(df) <- paste0(dimred, ncomponents) From 8efc0b75a160baaed7c7f2353f236e6bebff058d Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Mon, 3 Aug 2026 23:58:30 +0300 Subject: [PATCH 16/28] up --- R/plotOrdination.R | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 9f498d5f..f0f75a5b 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -368,13 +368,14 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), orig_attributes <- attributes(df) orig_attributes <- orig_attributes[ !names(orig_attributes) %in% c("dim", "dimnames") ] + # Subset the data to include specified columns to plot df <- df[, ncomponents] # Add colnames if they are not present if( is.null(colnames(df)) ){ colnames(df) <- paste0(dimred, ncomponents) } df <- df |> as.data.frame() - # Take only 2 specified columns + # Take names of the columns that will be plotted x_var <- colnames(df)[[1L]] y_var <- colnames(df)[[2L]] From 339cb1de20a3a9c1a64e5b796618436b6891f6b2 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 00:09:43 +0300 Subject: [PATCH 17/28] up --- R/plotOrdination.R | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index f0f75a5b..f142fcbe 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -392,8 +392,14 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), xlab <- x_var ylab <- y_var - eigen <- orig_attributes[["eig"]] - if( add.expl.var && !is.null(eigen) ){ + eigen_names <- c("eig", "percentVar") + if( add.expl.var && any(eigen_names %in% names(orig_attributes)) ){ + eigen_names <- eigen_names[ + eigen_names %in% names(orig_attributes)][[1L]] + eigen <- orig_attributes[[eigen_names]] + if( eigen_names %in% c("eig") ){ + eigen <- eigen * 100 + } xlab <- paste0( xlab, " (", round(eigen[ncomponents][[1L]], 1), From 680b548107d34c86e9644f8c79058c9abb3376e4 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 00:20:42 +0300 Subject: [PATCH 18/28] up --- R/plotOrdination.R | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index f142fcbe..db6509ac 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -751,6 +751,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), p, df, xlab = attributes(df)[["xlab"]], ylab = attributes(df)[["ylab"]], + coord.equal = TRUE, panel.by.eigen = TRUE, ...){ if( !.is_a_string(xlab) ){ @@ -759,6 +760,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !.is_a_string(ylab) ){ stop("'ylab' must be a single character value.", call. = FALSE) } + if( !.is_a_bool(coord.equal) ){ + stop("'coord.equal' must be TRUE or FALSE.", call. = FALSE) + } if( !.is_a_bool(panel.by.eigen) ){ stop("'panel.by.eigen' must be TRUE or FALSE.", call. = FALSE) } @@ -784,7 +788,9 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), # Use equal scaling on both axes so that one unit on the x-axis has the same # physical length as one unit on the y-axis. This preserves distances and # angles in the ordination and avoids visual distortion. - p <- p + coord_equal() + if( coord.equal ){ + p <- p + coord_equal() + } # Make the panel dimensions proportional to the eigenvalues so that the # relative lengths of the axes reflect the variation explained by each # ordination axis. From 21e324f96db44f1adb18533f42da5753b42f18a9 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 00:35:16 +0300 Subject: [PATCH 19/28] up --- R/plotOrdination.R | 32 +++++++++++++++++--------------- 1 file changed, 17 insertions(+), 15 deletions(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index db6509ac..dd884735 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -752,7 +752,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), xlab = attributes(df)[["xlab"]], ylab = attributes(df)[["ylab"]], coord.equal = TRUE, - panel.by.eigen = TRUE, + aspect.ratio = c("equal", "eigen", "free"), ...){ if( !.is_a_string(xlab) ){ stop("'xlab' must be a single character value.", call. = FALSE) @@ -763,9 +763,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !.is_a_bool(coord.equal) ){ stop("'coord.equal' must be TRUE or FALSE.", call. = FALSE) } - if( !.is_a_bool(panel.by.eigen) ){ - stop("'panel.by.eigen' must be TRUE or FALSE.", call. = FALSE) - } + match.arg(aspect.ratio) # p <- p + theme_classic() p <- p + labs(x = xlab, y = ylab) @@ -785,18 +783,22 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), p <- p + labs(linetype = attributes(df)[["linetype.by"]]) } - # Use equal scaling on both axes so that one unit on the x-axis has the same - # physical length as one unit on the y-axis. This preserves distances and - # angles in the ordination and avoids visual distortion. - if( coord.equal ){ + # Adjust axis length ratios + eigen_names <- c("eig", "percentVar") + if( aspect.ratio == "equal" ){ + # Use equal scaling on both axes so that one unit on the x-axis has the + # same physical length as one unit on the y-axis. This preserves + # distances and angles in the ordination and avoids visual distortion. p <- p + coord_equal() - } - # Make the panel dimensions proportional to the eigenvalues so that the - # relative lengths of the axes reflect the variation explained by each - # ordination axis. - if( panel.by.eigen && !is.null(attributes(df)[["eigen"]]) ){ - eig <- attributes(df)[["eigen"]] - p <- p + theme(aspect.ratio = eig[2] / eig[1]) + } else if ( aspect.ratio == "eigen" && + any(eigen_names %in% names(orig_attributes)) ){ + # Make the panel dimensions proportional to the eigenvalues so that the + # relative lengths of the axes reflect the variation explained by each + # ordination axis. + eigen_names <- eigen_names[ + eigen_names %in% names(orig_attributes)][[1L]] + eigen <- orig_attributes[[eigen_names]] + p <- p + theme(aspect.ratio = eigen[2] / eigen[1]) } # Adjust colors From 16d72ed5ddec82b9f4e654846b8d65b95cd588f5 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 14:54:18 +0300 Subject: [PATCH 20/28] up --- NEWS | 3 +++ 1 file changed, 3 insertions(+) diff --git a/NEWS b/NEWS index b55174bf..aa5b7d8a 100644 --- a/NEWS +++ b/NEWS @@ -51,3 +51,6 @@ Changes in version 1.17.x Changes in version 1.19.x + plotRDA: Now plotting works with interaction term (2025-11-09) + plotBoxplot: Added option to add p-values (2026-01-07) + +Changes in version 1.21.x ++ Added plotOrdination and plotJointRPCA (2026-08-04) From 471dacd28caf7c279bfa3913aafc0053497af2af Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 14:56:31 +0300 Subject: [PATCH 21/28] up --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 630b0a7d..85b69604 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -31,7 +31,7 @@ License: Artistic-2.0 | file LICENSE Encoding: UTF-8 LazyData: false Depends: - R (>= 4.0), + R (>= 4.1), ggplot2, ggraph (>= 2.0), mia (>= 1.13.0), From a0a019fa8ce23b29419380f484ac4d5e5c2f7b6e Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 14:59:03 +0300 Subject: [PATCH 22/28] fix --- R/plotOrdination.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index dd884735..035c36bb 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -763,7 +763,7 @@ setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), if( !.is_a_bool(coord.equal) ){ stop("'coord.equal' must be TRUE or FALSE.", call. = FALSE) } - match.arg(aspect.ratio) + aspect.ratio <- match.arg(aspect.ratio) # p <- p + theme_classic() p <- p + labs(x = xlab, y = ylab) From 3bc88e746911e5d77e4a3e89dcf13781c7fcebd9 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 15:26:11 +0300 Subject: [PATCH 23/28] up --- NAMESPACE | 3 +-- R/plotOrdination.R | 20 +++++++++++++++----- man/plotJointRPCA.Rd | 18 ++++++++++-------- man/plotOrdination.Rd | 17 +++++++++++++---- 4 files changed, 39 insertions(+), 19 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index dac933d2..6a80cf49 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -84,6 +84,7 @@ importFrom(SummarizedExperiment,colData) importFrom(SummarizedExperiment,rowData) importFrom(TreeSummarizedExperiment,colTree) importFrom(TreeSummarizedExperiment,colTreeNames) +importFrom(TreeSummarizedExperiment,rowLinks) importFrom(TreeSummarizedExperiment,rowTree) importFrom(TreeSummarizedExperiment,rowTreeNames) importFrom(ape,as.phylo) @@ -152,8 +153,6 @@ importFrom(ggplot2,scale_y_discrete) importFrom(ggplot2,theme) importFrom(ggplot2,theme_bw) importFrom(ggplot2,theme_classic) -importFrom(ggplot2,xlab) -importFrom(ggplot2,ylab) importFrom(ggrepel,geom_label_repel) importFrom(ggrepel,geom_text_repel) importFrom(ggtree,geom_cladelab) diff --git a/R/plotOrdination.R b/R/plotOrdination.R index 035c36bb..11ad2d0b 100644 --- a/R/plotOrdination.R +++ b/R/plotOrdination.R @@ -2,13 +2,23 @@ #' plotOrdination #' #' @title -#' Create ordination plot +#' Visualize ordination results #' #' @description -#' Ordinaton plotter +#' Creates a two-dimensional ordination plot from a reduced dimension result +#' stored in a +#' \code{\link[SingleCellExperiment:reducedDims]{SingleCellExperiment}} +#' object. Samples can be coloured, filled, shaped, sized, grouped or faceted +#' using sample metadata or feature abundances, and additional graphical +#' elements such as confidence ellipses, centroids, vectors and density +#' estimates can be added. #' #' @details -#' Creates ordination plot +#' This function provides a unified interface for visualizing ordination methods +#' such as PCA, PCoA, MDS, t-SNE, UMAP, RDA and CCA. The plotted coordinates are +#' retrieved from a reduced dimension result stored in +#' \code{reducedDim(x)}. +#' #' #' @return #' A \code{ggplot2} object. @@ -249,8 +259,8 @@ NULL #' @rdname plotOrdination #' @export setMethod("plotOrdination", signature = c(x = "SingleCellExperiment"), - function(x, dimred, colour.by = color.by, color.by = NULL, ...){ - args <- .check_ordination_input(x, dimred, colour.by = colour.by, ...) + function(x, dimred, ...){ + args <- .check_ordination_input(x, dimred, ...) df <- do.call(.get_ordination_data, args) p <- .ordination_plotter(df, ...) return(p) diff --git a/man/plotJointRPCA.Rd b/man/plotJointRPCA.Rd index f0adbc7a..9162ebf3 100644 --- a/man/plotJointRPCA.Rd +++ b/man/plotJointRPCA.Rd @@ -10,7 +10,7 @@ plotJointRPCA(x, ...) \S4method{plotJointRPCA}{MultiAssayExperiment}(x, dimred, ...) -\S4method{plotJointRPCA}{SingleCellExperiment}(x, dimred, ...) +\S4method{plotJointRPCA}{SingleCellExperiment}(x, dimred, add.vectors = TRUE, ...) } \arguments{ \item{x}{A @@ -24,16 +24,18 @@ object containing Joint-RPCA results.} plotting. Commonly used Joint-RPCA-specific arguments include: \itemize{ -\item \code{add.vectors}: \code{logical(1)} or \code{character}. If -\code{TRUE}, feature loading vectors are added. Alternatively, a character -vector can be supplied to display only features whose names match the given -pattern(s). (Default: \code{TRUE}) \item \code{ntop}: \code{integer(1)} or \code{NULL}. Maximum number of loading vectors shown per data layer. The vectors with the largest lengths are retained. Set to \code{NULL} to display all vectors. (Default: \code{10}) }} -\item{dimred}{\code{character(1)} specifying the name of the Joint-RPCA result to plot.} +\item{dimred}{\code{character(1)}: Specifies the name of the Joint-RPCA +result to plot.} + +\item{add.vectors:}{\code{logical(1)} or \code{character}. If +\code{TRUE}, feature loading vectors are added. Alternatively, a character +vector can be supplied to display only features whose names match the given +pattern(s). (Default: \code{TRUE})} } \value{ A \code{ggplot2} object. @@ -98,14 +100,14 @@ plotJointRPCA( ntop = NULL ) -## Display only selected feature vectors +# Display only selected feature vectors plotJointRPCA( mae, "JointRPCA", add.vectors = c("Bacteroides", "Roseburia") ) -## Use boxed labels without repelling +# Use boxed labels without repelling plotJointRPCA( mae, "JointRPCA", diff --git a/man/plotOrdination.Rd b/man/plotOrdination.Rd index 7be443e9..c005b8a9 100644 --- a/man/plotOrdination.Rd +++ b/man/plotOrdination.Rd @@ -3,11 +3,11 @@ \name{plotOrdination} \alias{plotOrdination} \alias{plotOrdination,SingleCellExperiment-method} -\title{Create ordination plot} +\title{Visualize ordination results} \usage{ plotOrdination(x, ...) -\S4method{plotOrdination}{SingleCellExperiment}(x, dimred, colour.by = color.by, color.by = NULL, ...) +\S4method{plotOrdination}{SingleCellExperiment}(x, dimred, ...) } \arguments{ \item{x}{a @@ -138,10 +138,19 @@ stored in \code{reducedDim(x)} to visualize.} A \code{ggplot2} object. } \description{ -Ordinaton plotter +Creates a two-dimensional ordination plot from a reduced dimension result +stored in a +\code{\link[SingleCellExperiment:reducedDims]{SingleCellExperiment}} +object. Samples can be coloured, filled, shaped, sized, grouped or faceted +using sample metadata or feature abundances, and additional graphical +elements such as confidence ellipses, centroids, vectors and density +estimates can be added. } \details{ -Creates ordination plot +This function provides a unified interface for visualizing ordination methods +such as PCA, PCoA, MDS, t-SNE, UMAP, RDA and CCA. The plotted coordinates are +retrieved from a reduced dimension result stored in +\code{reducedDim(x)}. } \examples{ data("Tito2024QMP") From cca18ff97a9c25a208e7c7c9531ab6a0dabbbdfd Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 15:32:03 +0300 Subject: [PATCH 24/28] up --- R/plotCCA.R | 13 ------------- 1 file changed, 13 deletions(-) diff --git a/R/plotCCA.R b/R/plotCCA.R index fdceaeec..80c851ed 100644 --- a/R/plotCCA.R +++ b/R/plotCCA.R @@ -193,19 +193,6 @@ setMethod("plotRDA", signature = c(x = "SingleCellExperiment"), ################################ HELP FUNCTIONS ################################ -# Construct TreeSE from matrix to pass it to downstream functions. It is useful -# for instance if get* functios was used instead of add*. -#' @importFrom S4Vectors SimpleList -.rda2tse <- function(object) { - # Convert rda/cca object to TreeSE - object <- TreeSummarizedExperiment( - assays = SimpleList(counts = matrix( - ncol = nrow(object), dimnames = list(NULL, rownames(object)))), - reducedDims = list(RDA = object) - ) - return(object) -} - # This function retrieves data for creating vectors. Moreover, it wrangles the # vector data and controls what information is added to vector text or labels. .get_rda_vector_data <- function( From ed30c214058143e00c21157bea77edf7d5719750 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 15:33:54 +0300 Subject: [PATCH 25/28] up --- R/plotJointRPCA.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/plotJointRPCA.R b/R/plotJointRPCA.R index 38aa222b..2ccd417f 100644 --- a/R/plotJointRPCA.R +++ b/R/plotJointRPCA.R @@ -118,7 +118,7 @@ NULL #' @export setMethod("plotJointRPCA", signature = c(x = "MultiAssayExperiment"), function(x, dimred, ...){ - mia:::.check_metadata_present(dimred, x) + .check_metadata_present(dimred, x) # Construct TreeSE from results so that we can use TreeSE function # (and thus plotOrdination) res <- metadata(x)[[dimred]] From 571ee6a77e98e043bd0a8ebbd528df5128076473 Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 15:35:01 +0300 Subject: [PATCH 26/28] up --- R/plotJointRPCA.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/plotJointRPCA.R b/R/plotJointRPCA.R index 2ccd417f..232d778f 100644 --- a/R/plotJointRPCA.R +++ b/R/plotJointRPCA.R @@ -167,8 +167,7 @@ setMethod("plotJointRPCA", signature = c(x = "SingleCellExperiment"), # Get feature loadings that will be plotted as vectors #' @importFrom dplyr group_by slice_max ungroup .get_joint_rpca_vector_data <- function( - tse, reduced_dim, add.vectors, ignore.case = FALSE, ntop = 10, - ...){ + tse, reduced_dim, add.vectors, ignore.case = FALSE, ntop = 10, ...){ if( !.is_a_bool(ignore.case) ){ stop("'ignore.case' must be TRUE or FALSE.", call. = FALSE) } From fef1a4f21ca5a5a517a1fd8c5838e47e85c2106a Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Tue, 4 Aug 2026 16:52:10 +0300 Subject: [PATCH 27/28] Fix tests. Default choices have changed for plotRDA --- tests/testthat/test-plotCCA.R | 25 +++++++++++++------------ 1 file changed, 13 insertions(+), 12 deletions(-) diff --git a/tests/testthat/test-plotCCA.R b/tests/testthat/test-plotCCA.R index 798278ad..53393bdd 100644 --- a/tests/testthat/test-plotCCA.R +++ b/tests/testthat/test-plotCCA.R @@ -10,16 +10,15 @@ test_that("plot RDA/CCA", { # Run/calculate RDA tse <- addRDA(tse, assay.type = "counts", formula = assay ~ patient_status + cohort) - rda <- getRDA(tse, assay.type = "counts", formula = assay ~ patient_status + cohort) # Minimal functionality expect_no_error(plotRDA(tse, "RDA")) # Wrong-entry scenarios - expect_error(plotRDA(tse, "RDA", colour_by = "wrong colname")) - expect_error(plotRDA(tse, "RDA", colour_by = "cohort", shape_by = "wrong colname")) + expect_error(plotRDA(tse, "RDA", colour.by = "wrong colname")) + expect_error(plotRDA(tse, "RDA", colour.by = "cohort", shape.by = "wrong colname")) expect_error(plotRDA(tse, "RDA", add.ellipse = "invalid value"), - "'add.ellipse' must be one of c(TRUE, FALSE, 'fill', 'color').", + "'add.ellipse' must be TRUE or FALSE.", fixed = TRUE) expect_error(plotRDA(tse, "RDA", add.significance = "invalid value"), "'add.significance' must be TRUE or FALSE.") @@ -34,12 +33,12 @@ test_that("plot RDA/CCA", { ### 2). TEST plot layers ### - el_true <- plotRDA(tse, "RDA", colour_by = "patient_status") - el_false <- plotRDA(tse, "RDA", colour_by = "patient_status", add.ellipse = FALSE) - el_col <- plotRDA(tse, "RDA", colour_by = "patient_status", add.ellipse = "colour") - el_fill <- plotRDA(tse, "RDA", colour_by = "patient_status", add.ellipse = "fill") + el_true <- plotRDA(tse, "RDA", colour.by = "patient_status", add.ellipse = TRUE) + el_false <- plotRDA(tse, "RDA", colour.by = "patient_status", add.ellipse = FALSE) + el_col <- plotRDA(tse, "RDA", colour.by = "patient_status", add.ellipse = TRUE) + el_fill <- plotRDA(tse, "RDA", fill.by = "patient_status", add.ellipse = TRUE) expect_warning( - vec_false <- plotRDA(tse, "RDA", colour_by = "patient_status", add.vectors = FALSE) + vec_false <- plotRDA(tse, "RDA", colour.by = "patient_status", add.ellipse = TRUE, add.vectors = FALSE) ) # Filled ellipse has one more layer than no ellipse plot expect_equal(length(ggplot_build(el_true)[["data"]]), 4) @@ -51,7 +50,9 @@ test_that("plot RDA/CCA", { expect_false(all(ggplot_build(el_fill)[["data"]][[2]][["alpha"]] == 0)) # Check ggplot aesthetics - p_aes <- plotRDA(tse, "RDA", colour_by = "patient_status", ellipse.alpha = 0.5, + p_aes <- plotRDA(tse, "RDA", + colour.by = "patient_status", fill.by = "patient_status", + add.ellipse = TRUE, ellipse.alpha = 0.5, ellipse.linewidth = 0.2, ellipse.linetype = 3, vec.size = 0.6, vec.colour = "red", vec.linetype = 2, arrow.size = 0.15, label.colour = "blue", label.size = 5) @@ -72,8 +73,8 @@ test_that("plot RDA/CCA", { # expect_true(arrow_size == 0.15) # Vector or label text - p_vec <- plotRDA(tse, "RDA", colour_by = "patient_status", vec.text = TRUE) - p_lab <- plotRDA(tse, "RDA", colour_by = "patient_status", vec.text = FALSE) + p_vec <- plotRDA(tse, "RDA", colour.by = "patient_status", fill.by = "patient_status", add.ellipse = TRUE, vec.text = TRUE) + p_lab <- plotRDA(tse, "RDA", colour.by = "patient_status", fill.by = "patient_status", add.ellipse = TRUE, vec.text = FALSE) # There must be label column in both expect_true( "label" %in% (ggplot_build(p_vec)[["data"]][[4]] |> names()) ) expect_true( "label" %in% (ggplot_build(p_lab)[["data"]][[4]] |> names()) ) From 8f79a503b4fa235806d331c5848cb20b3ee0a66f Mon Sep 17 00:00:00 2001 From: Tuomas Borman Date: Wed, 5 Aug 2026 11:36:59 +0300 Subject: [PATCH 28/28] up --- NAMESPACE | 1 - R/plotJointRPCA.R | 2 +- man/plotJointRPCA.Rd | 2 +- 3 files changed, 2 insertions(+), 3 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 6a80cf49..fa0a1614 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -72,7 +72,6 @@ importFrom(BiocParallel,bpstop) importFrom(DelayedArray,rowMeans) importFrom(DelayedArray,rowSums) importFrom(DirichletMultinomial,mixture) -importFrom(S4Vectors,SimpleList) importFrom(S4Vectors,metadata) importFrom(S4Vectors,unfactor) importFrom(SingleCellExperiment,reducedDim) diff --git a/R/plotJointRPCA.R b/R/plotJointRPCA.R index 232d778f..3611bd19 100644 --- a/R/plotJointRPCA.R +++ b/R/plotJointRPCA.R @@ -33,7 +33,7 @@ #' @param dimred \code{character(1)}: Specifies the name of the Joint-RPCA #' result to plot. #' -#' @param add.vectors: \code{logical(1)} or \code{character}. If +#' @param add.vectors \code{logical(1)} or \code{character}. If #' \code{TRUE}, feature loading vectors are added. Alternatively, a character #' vector can be supplied to display only features whose names match the given #' pattern(s). (Default: \code{TRUE}) diff --git a/man/plotJointRPCA.Rd b/man/plotJointRPCA.Rd index 9162ebf3..7a751124 100644 --- a/man/plotJointRPCA.Rd +++ b/man/plotJointRPCA.Rd @@ -32,7 +32,7 @@ are retained. Set to \code{NULL} to display all vectors. (Default: \code{10}) \item{dimred}{\code{character(1)}: Specifies the name of the Joint-RPCA result to plot.} -\item{add.vectors:}{\code{logical(1)} or \code{character}. If +\item{add.vectors}{\code{logical(1)} or \code{character}. If \code{TRUE}, feature loading vectors are added. Alternatively, a character vector can be supplied to display only features whose names match the given pattern(s). (Default: \code{TRUE})}