From 4bed6bf7d82e69f3a65d3db84fabec0cfbf44b49 Mon Sep 17 00:00:00 2001 From: phauchamps Date: Thu, 4 Jan 2024 15:10:51 +0100 Subject: [PATCH] getChannelSummaryStats() : limit case with one flowFrame was still not working according to expectations :-) --- R/stats.R | 29 ++++++++++++++++------------- man/getChannelSummaryStats.Rd | 5 ++++- tests/testthat/test-stats.R | 33 +++++++++++++++++++-------------- 3 files changed, 39 insertions(+), 28 deletions(-) diff --git a/R/stats.R b/R/stats.R index 0d3e53d..7854341 100644 --- a/R/stats.R +++ b/R/stats.R @@ -1009,8 +1009,11 @@ pairwiseEMDDist <- function( #' This list can be named, in that case, these names will be transfered to the #' returned value. #' @param ... additional parameters passed to `getEMDDist()` -#' @return a matrix of which the columns are the channel statistics -#' for all flowFrames of the flowSet. +#' @return a list of named statistic matrices. +#' In each stat matrix, the columns are the channel statistics +#' for all flowFrames of the flowSet. +#' Exception: if only one stat function (and not a list) is passed in +#' `statFUNs`, the return value is simplified to the stat matrix itself. #' @importFrom CytoPipeline areSignalCols #' @export #' @@ -1132,20 +1135,20 @@ getChannelSummaryStats <- function( chRes <- flowCore::fsApply( fs, FUN = function(ff){ - statFUNList[[fu]](flowCore::exprs(ff)[, ch], - na.rm = TRUE) + statFUNList[[fu]]( + flowCore::exprs(ff)[, ch, drop = FALSE], + na.rm = TRUE) }) }, FUN.VALUE = numeric(length = nFF)) - if (is.vector(statList[[fu]])) { - names(statList[[fu]]) <- channelNames - } else { # matrix - colnames(statList[[fu]]) <- channelNames - rowNames <- flowCore::pData(fs)$name - if (!is.null(rowNames)) { - rownames(statList[[fu]]) <- rowNames - } + if (!is.matrix(statList[[fu]]) ) { + statList[[fu]] <- matrix(statList[[fu]], nrow = 1) + } + colnames(statList[[fu]]) <- channelNames + rowNames <- flowCore::pData(fs)$name + if (!is.null(rowNames)) { + rownames(statList[[fu]]) <- rowNames } } @@ -1153,7 +1156,7 @@ getChannelSummaryStats <- function( names(statList) <- names(statFUNList) # if only one stat function => unlist to return one single matrix - if (nStats == 1){ + if (nStats == 1 && !is.list(statFUNs)){ statList <- statList[[1]] } diff --git a/man/getChannelSummaryStats.Rd b/man/getChannelSummaryStats.Rd index 6b4e1f9..aa12020 100644 --- a/man/getChannelSummaryStats.Rd +++ b/man/getChannelSummaryStats.Rd @@ -37,8 +37,11 @@ after each single distance calculation} \item{...}{additional parameters passed to \code{getEMDDist()}} } \value{ -a matrix of which the columns are the channel statistics +a list of named statistic matrices. +In each stat matrix, the columns are the channel statistics for all flowFrames of the flowSet. +Exception: if only one stat function (and not a list) is passed in +\code{statFUNs}, the return value is simplified to the stat matrix itself. } \description{ Calculate a summary statistic of some channels of diff --git a/tests/testthat/test-stats.R b/tests/testthat/test-stats.R index 6354873..ea0042c 100644 --- a/tests/testthat/test-stats.R +++ b/tests/testthat/test-stats.R @@ -716,6 +716,8 @@ test_that("getChannelSummaryStats works", { statFUNs = stats::median, verbose = FALSE) + + expect_equal(unname(rownames(ret)), flowCore::sampleNames(fsAll)) expect_equal(unname(colnames(ret)), channelsOrMarkers) expect_equal(unname(ret[2,1]), 1.941572819633120) @@ -729,20 +731,19 @@ test_that("getChannelSummaryStats works", { verbose = FALSE) expect_equal(names(ret), c("mean", "std.dev")) + expect_equal(unname(rownames(ret[[1]])), flowCore::sampleNames(fsAll)) expect_equal(unname(colnames(ret[[1]])), channelsOrMarkers) expect_equal(unname(ret[[1]][2,1]), 1.961771149533659) expect_equal(unname(ret[[1]][3,2]), 1.393034241892733) expect_equal(unname(ret[[1]][4,3]), 1.907561024702794) - + expect_equal(unname(rownames(ret[[2]])), flowCore::sampleNames(fsAll)) expect_equal(unname(colnames(ret[[2]])), channelsOrMarkers) - expect_equal(unname(ret[[2]][2,1]), 0.5942291770870205) expect_equal(unname(ret[[2]][3,2]), 0.7352746696905406) expect_equal(unname(ret[[2]][4,3]), 0.8420740225208073) - # test with only one flow frame => should return a vector - # instead of a matrix + # test with only one flow frame ret <- getChannelSummaryStats( fsAll[1], channels = channelsOrMarkers, @@ -750,8 +751,10 @@ test_that("getChannelSummaryStats works", { verbose = FALSE) expect_equal(names(ret), c("mean", "std.dev")) - expect_equal(unname(names(ret[[1]])), channelsOrMarkers) - expect_equal(unname(names(ret[[2]])), channelsOrMarkers) + expect_equal(unname(rownames(ret[[1]])), "Donor1") + expect_equal(unname(rownames(ret[[2]])), "Donor1") + expect_equal(unname(colnames(ret[[1]])), channelsOrMarkers) + expect_equal(unname(colnames(ret[[2]])), channelsOrMarkers) expect_equal(unname(ret[[1]][1]), 1.900298) expect_equal(unname(ret[[1]][2]), 1.39186533) expect_equal(unname(ret[[1]][3]), 1.8544648) @@ -766,10 +769,11 @@ test_that("getChannelSummaryStats works", { statFUNs = list("mean" = mean), verbose = FALSE) - expect_equal(names(ret), channelsOrMarkers) - expect_equal(unname(ret[1]), 1.900298) - expect_equal(unname(ret[2]), 1.39186533) - expect_equal(unname(ret[3]), 1.8544648) + expect_equal(as.character(rownames(ret[[1]])), "Donor1") + expect_equal(unname(colnames(ret[[1]])), channelsOrMarkers) + expect_equal(unname(ret[[1]][1,1]), 1.900298) + expect_equal(unname(ret[[1]][1,2]), 1.39186533) + expect_equal(unname(ret[[1]][1,3]), 1.8544648) # one flow frame, one stat (bis with direct impact of stat FUN) ret <- getChannelSummaryStats( @@ -778,10 +782,11 @@ test_that("getChannelSummaryStats works", { statFUNs = mean, verbose = FALSE) - expect_equal(names(ret), channelsOrMarkers) - expect_equal(unname(ret[1]), 1.900298) - expect_equal(unname(ret[2]), 1.39186533) - expect_equal(unname(ret[3]), 1.8544648) + expect_equal(as.character(rownames(ret)), "Donor1") + expect_equal(unname(colnames(ret)), channelsOrMarkers) + expect_equal(unname(ret[1,1]), 1.900298) + expect_equal(unname(ret[1,2]), 1.39186533) + expect_equal(unname(ret[1,3]), 1.8544648) }) test_that("computeMetricMDS works", {