From 4a66f6c2b3c481319105dccfdc04ffb33d64c9e4 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Fri, 17 Oct 2025 23:57:41 +0200 Subject: [PATCH 1/7] Add full SummarizedExperiment support --- DESCRIPTION | 3 +- R/bbdml.R | 8 + R/contrastsTest.R | 10 + R/convert_sumexp.R | 22 ++ R/differentialTest.R | 15 +- R/plot_differentialTest.R | 16 + R/warn_phyloseq.R | 7 + tests/testthat/test_testsRequiringSE.R | 506 +++++++++++++++++++++++++ 8 files changed, 585 insertions(+), 2 deletions(-) create mode 100644 R/convert_sumexp.R create mode 100644 tests/testthat/test_testsRequiringSE.R diff --git a/DESCRIPTION b/DESCRIPTION index 60dd9fa..132c868 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -26,5 +26,6 @@ Suggests: knitr, slam, R.rsp, optimx, - phyloseq + phyloseq, + SummarizedExperiment VignetteBuilder: knitr diff --git a/R/bbdml.R b/R/bbdml.R index 1032671..df1089c 100644 --- a/R/bbdml.R +++ b/R/bbdml.R @@ -60,6 +60,14 @@ bbdml <- function(formula, phi.formula, data, formula <- stats::update(formula, cbind(W, M - W) ~ .) } + # Convert SummarizedExperiment objects + if (inherits(data, "SummarizedExperiment")) { + selection <- all.vars(formula)[1] + data <- convert_sumexp(data, select = selection) + # Update formula to match convert_phylo specification + formula <- stats::update(formula, cbind(W, M - W) ~ .) + } + # Record call call <- match.call(expand.dots = FALSE) # Record mu link diff --git a/R/contrastsTest.R b/R/contrastsTest.R index 396cd86..13440e4 100644 --- a/R/contrastsTest.R +++ b/R/contrastsTest.R @@ -95,6 +95,14 @@ contrastsTest <- function(formula, phi.formula, } else { warn_phyloseq() } + } else if (inherits(data, "SummarizedExperiment")) { + if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + # Set up response + taxanames <- row.names(SummarizedExperiment::rowData(data)) + sample_data <- SummarizedExperiment::colData(data) + } else { + warn_sumexp() + } } else if (is.matrix(data) || is.data.frame(data)) { # # use phyloseq @@ -147,6 +155,8 @@ contrastsTest <- function(formula, phi.formula, # Subset data to only select that taxa if ("phyloseq" %in% class(data)) { data_i <- convert_phylo(data, select = taxanames[i]) + } else if (inherits(data, "SummarizedExperiment")) { + data_i <- convert_sumexp(data, select = taxanames[i]) } else { response_i <- data.frame(W = data[, taxanames[i]], M = M) data_i <- cbind(response_i, sample_data) diff --git a/R/convert_sumexp.R b/R/convert_sumexp.R new file mode 100644 index 0000000..48e2822 --- /dev/null +++ b/R/convert_sumexp.R @@ -0,0 +1,22 @@ +#' Function to subset and convert SummarizedExperiment data +#' +#' @param data a \code{SummarizedExperiment} object +#' @param select Name of OTU or taxa to select, must match taxa name in \code{data} +#' +#' @return A \code{data.frame} object, with elements \code{W} as the observed counts, \code{M} as the sequencing depth, and the sample data with their original names. +#' +#' @export +convert_sumexp <- function(data, select) { + if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + subsamp <- data[select, ] + W_tmp <- matrix(t(SummarizedExperiment::assay(subsamp)), ncol = 1) + + + out <- data.frame(W = W_tmp, + M = colSums(SummarizedExperiment::assay(data)), + SummarizedExperiment::colData(subsamp)) + return(out) + } else { + warn_sumexp() + } +} diff --git a/R/differentialTest.R b/R/differentialTest.R index 3843bd7..bef29a5 100644 --- a/R/differentialTest.R +++ b/R/differentialTest.R @@ -95,6 +95,15 @@ differentialTest <- function(formula, phi.formula, } else { warn_phyloseq() } + } else if (inherits(data, "SummarizedExperiment")) { + if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + # Set up response + taxanames <- row.names(SummarizedExperiment::rowData(data)) + sample_data <- data.frame(SummarizedExperiment::colData(data)) + } else { + warn_sumexp() + } + } else if (is.matrix(data) || is.data.frame(data)) { # # use phyloseq @@ -118,7 +127,7 @@ differentialTest <- function(formula, phi.formula, M <- rowSums(data) } else { - stop("Input must be either data frame, matrix, or phyloseq object!") + stop("Input must be either data frame, matrix, phyloseq object or SummarizedExperiment!") } # Set up output @@ -158,6 +167,8 @@ differentialTest <- function(formula, phi.formula, # Subset data to only select that taxa if ("phyloseq" %in% class(data)) { data_i <- convert_phylo(data, select = taxanames[i]) + } else if (inherits(data, "SummarizedExperiment")) { + data_i <- convert_sumexp(data, select = taxanames[i]) } else { response_i <- data.frame(W = data[, taxanames[i]], M = M) data_i <- cbind(response_i, sample_data) @@ -275,6 +286,8 @@ We *strongly recommend* running `bbdml` on a single taxon (especially before pos i <- (try_only[!(try_only %in% ind_disc)])[1] if ("phyloseq" %in% class(data)) { data_i <- convert_phylo(data, select = taxanames[i]) + } else if (inherits(data, "SummarizedExperiment")) { + data_i <- convert_sumexp(data, select = taxanames[i]) } else { response_i <- data.frame(W = data[, taxanames[i]], M = M) data_i <- cbind(response_i, sample_data) diff --git a/R/plot_differentialTest.R b/R/plot_differentialTest.R index 701f625..1553ccc 100644 --- a/R/plot_differentialTest.R +++ b/R/plot_differentialTest.R @@ -35,6 +35,22 @@ plot.differentialTest <- function(x, level = NULL, data_only = FALSE, ...) { } } } + + if (inherits(x$data, "SummarizedExperiment")) { + if (!(nrow(SummarizedExperiment::rowData(x$data)) == 0)) { + if (is.null(level)) { + signif_taxa = apply(SummarizedExperiment::rowData(x$data)[signif_taxa,], 1, function(x) {paste(stats::na.omit(x), collapse = '_')}) + } else { + signif_taxa = sapply(SummarizedExperiment::rowData(x$data)[signif_taxa, level], function(x) {paste(stats::na.omit(x), collapse = '_')}) + } + + if (length(unique(signif_taxa)) != length(unique(x$significant_taxa))) { + # Make sure if repeated taxa add unique otu identifiers + signif_taxa <- paste0(signif_taxa, " (", x$significant_taxa, ")") + } + } + } + if (length(x$significant_models) != 0) { var_per_mod <- length(x$restrictions_DA) + length(x$restrictions_DV) total_var_count <- length(signif_taxa) * var_per_mod diff --git a/R/warn_phyloseq.R b/R/warn_phyloseq.R index c8c9d1b..7757722 100644 --- a/R/warn_phyloseq.R +++ b/R/warn_phyloseq.R @@ -4,3 +4,10 @@ warn_phyloseq <- function() { stop("You are trying to use a `phyloseq` data object or `phyloseq` helper function without having the `phyloseq` package installed. Please either install the package or use a standard data frame.") } + +#' Function to throw error if the `SummarizedExperiment` package is called but it is not installed +#' +#' @export +warn_sumexp <- function() { + stop("You are trying to use a `SummarizedExperiment` data object or `SummarizedExperiment` helper function without having the `SummarizedExperiment` package installed. Please either install the package or use a standard data frame.") +} diff --git a/tests/testthat/test_testsRequiringSE.R b/tests/testthat/test_testsRequiringSE.R new file mode 100644 index 0000000..91a0882 --- /dev/null +++ b/tests/testthat/test_testsRequiringSE.R @@ -0,0 +1,506 @@ +library(corncob) +context("Tests that use SummarizedExperiment") + +data("soil_phylo_sample") +data("soil_phylo_otu") +data("soil_phylo_taxa") +data("soil_phylum_small_sample") +data("soil_phylum_small_otu") +data("soil_phylum_contrasts_sample") +data("soil_phylum_contrasts_otu") + +if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + soil_phylo_ps <- phyloseq::phyloseq( + phyloseq::sample_data(soil_phylo_sample), + phyloseq::otu_table(soil_phylo_otu, taxa_are_rows = TRUE), + phyloseq::tax_table(soil_phylo_taxa) + ) + + # make SummarizedExperiment objects from data + soil_phylo <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil_phylo_sample)), + assays = list(counts = data.frame(phyloseq::otu_table(soil_phylo_otu, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(soil_phylo_taxa)) + ) + soil_phylum_small <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil_phylum_small_sample)), + assays = list(counts = data.frame(phyloseq::otu_table(soil_phylum_small_otu, taxa_are_rows = TRUE))) + ) + soil_phylum_small_contrasts <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil_phylum_contrasts_sample)), + assays = list(counts = data.frame(phyloseq::otu_table(soil_phylum_contrasts_otu, taxa_are_rows = TRUE))) + ) + + # tests for bbdml() + test_that("bbdml works with SummarizedExperiment object", { + out_phylo <- bbdml(formula = Proteobacteria ~ 1, + phi.formula = ~ 1, + data = soil_phylum_small) + expect_is(out_phylo, "bbdml") + }) + + # tests for contrastsTest() + set.seed(1) + limma_install <- try(find.package("limma"), silent = TRUE) + if (!(inherits(limma_install, "try-error"))) { + temp <- contrastsTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + contrasts_DA = list("DayAmdmt21 - DayAmdmt11", + "DayAmdmt22 - DayAmdmt21"), + data = soil_phylum_small_contrasts, + fdr_cutoff = 0.05) + } + + test_that("contrastTest works", { + if (!(inherits(limma_install, "try-error"))) { + expect_is(temp, "contrastsTest") + } else { + expect_error(contrastsTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + contrasts_DA = list("DayAmdmt21 - DayAmdmt11", + "DayAmdmt22 - DayAmdmt21"), + data = soil_phylum_small_contrasts, + fdr_cutoff = 0.05), + "If you would like to test contrasts, please install the `limma` package, available through Bioconductor.") + } + }) + + # tests for getRestrictionTerms() + test_that("getRestrictionTerms works for phyloseq object", { + soil <- soil_phylo_ps %>% + phyloseq::subset_samples(DayAmdmt %in% c(11,21)) %>% + phyloseq::tax_glom("Phylum") + + soil <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil)), + assays = list(counts = data.frame(phyloseq::otu_table(soil, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(soil)) + ) + + mod1 <- bbdml(formula = OTU.1 ~ Day*Plants, + phi.formula = ~ Plants, + data = soil) + + mod2 <- bbdml(formula = OTU.1 ~Day - 1, + phi.formula = ~ Plants - 1, + data = soil) + + tmp <- corncob:::getRestrictionTerms(mod1,mod2) + expect_equal(tmp$mu, c(1,3,4)) + expect_true(tmp$phi == 5) + }) + + # tests for differentialTest() + set.seed(1) + soil <- phyloseq::subset_samples(soil_phylo_ps, DayAmdmt %in% c(11,21)) + subsoil <- phyloseq::prune_taxa(x = soil, taxa = rownames(phyloseq::otu_table(soil))[301:325]) + subsoil_ps <- phyloseq::prune_taxa(x = soil, taxa = rownames(phyloseq::otu_table(soil))[301:325]) + # including 7027 for all zeros + subsoil_disc <- phyloseq::prune_taxa(x = soil, taxa = rownames(phyloseq::otu_table(soil))[c(3001:3024, 7027)]) + + subsoil_nonint <- subsoil + phyloseq::otu_table(subsoil_nonint) <- phyloseq::otu_table(subsoil_nonint) + rexp(1) + + soil <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil)), + assays = list(counts = data.frame(phyloseq::otu_table(soil, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(soil)) + ) + + subsoil <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(subsoil)), + assays = list(counts = data.frame(phyloseq::otu_table(subsoil, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(subsoil)) + ) + + subsoil_disc <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(subsoil_disc)), + assays = list(counts = data.frame(phyloseq::otu_table(subsoil_disc, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(subsoil_disc)) + ) + + subsoil_nonint <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(subsoil_nonint)), + assays = list(counts = data.frame(phyloseq::otu_table(subsoil_nonint, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(subsoil_nonint)) + ) + + temp <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp2 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + filter_discriminant = FALSE, + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp3 <- suppressWarnings( + differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_disc, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + # Record of expected warnings below # + + # Warning in print.bbdml(mod) : + # This model is based on a discriminant taxa. + # You may see NAs in the model summary because Wald testing is invalid. + # Likelihood ratio testing can be used, but valid standard errors cannot be calculated. + + # Warning in waldt(object) : + # Singular Hessian! Cannot calculate p-values in this setting. + ) + + temp_badinits1 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(Inf, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp_badinits2 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(Inf, 2))) + + # Add this to cause some warnings and check those + temp_noinit_sing <- differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, boot = FALSE, test = "LRT", + data = subsoil) + + temp_noinit <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, boot = FALSE, test = "LRT", + data = subsoil) + + temp_sing <- differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 4))) + + temp_sing_ps <- differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_ps, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 4))) + + + temp_badinits3 <- differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(Inf, 4))) + + mydat <- SummarizedExperiment::assay(subsoil) + mysampdat <- SummarizedExperiment::colData(subsoil) + + temp_nonphylo <- differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = mydat, boot = FALSE, test = "LRT", + sample_data = mysampdat, + inits = rbind(rep(.01, 4))) + + temp_wald <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "Wald", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp_pblrt <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = TRUE, B = 5, test = "LRT", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp_pbwald <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = TRUE, B = 5, test = "Wald", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2))) + + temp_pbwald_robust <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = TRUE, B = 5, test = "Wald", + inits = rbind(rep(.01, 6)), + inits_null = rbind(rep(0.01, 2)), + robust = TRUE) + + + test_that("differentialTest works", { + expect_is(temp, "differentialTest") + expect_is(temp2, "differentialTest") + expect_is(temp3, "differentialTest") + expect_is(temp_wald, "differentialTest") + expect_is(temp_pbwald, "differentialTest") + expect_is(temp_pblrt, "differentialTest") + expect_is(temp, "differentialTest") + expect_is(temp_sing, "differentialTest") + expect_is(temp_nonphylo, "differentialTest") + expect_is(temp_noinit, "differentialTest") + expect_is(temp_noinit_sing, "differentialTest") + expect_is(temp_badinits1, "differentialTest") + expect_is(temp_badinits2, "differentialTest") + expect_is(temp_badinits3, "differentialTest") + expect_is(temp_pbwald_robust, "differentialTest") + }) + + + test_that("differentialTest S3 methods", { + expect_is(plot(temp), "ggplot") + expect_is(plot(temp, level = c("Order", "Class")), "ggplot") + expect_is(plot(temp, level = "Kingdom"), "ggplot") + expect_output(expect_null(print(temp))) + }) + + test_that("differentialTest works without SummarizedExperiment", { + expect_true(all.equal(temp_sing$p, temp_nonphylo$p)) + }) + + test_that("differentialTest phyloseq and SummarizedExperiment are identical", { + expect_true(all.equal(temp_sing$p, temp_sing_ps$p)) + }) + + + test_that("requires data frame, matrix, phyloseq, or SummarizedExperiment", { + expect_error(differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = c(1,2,3), boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 4)))) + }) + + test_that("inits require correct length", { + expect_error(differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits = rbind(rep(.01, 6)))) + expect_error(differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + inits_null = rbind(rep(.01, 4)))) + }) + + test_that("try_only works", { + expect_is(differentialTest(formula = ~ DayAmdmt, + phi.formula = ~ DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + try_only = 1:2), "differentialTest") + }) + + test_that("overspecification error message", { + expect_error(differentialTest(formula = ~ Plants*Day*Amdmt, + phi.formula = ~ Plants*Day*Amdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, boot = FALSE, test = "LRT", + try_only = 1:2)) + }) + + test_that("differentialTest does NAs correctly", { + expect_equal(length(temp$p), 25) + expect_equal(length(temp$p_fdr), 25) + expect_equal(length(temp3$p), 25) + expect_equal(length(temp3$p_fdr), 25) + }) + + + + test_that("differentialTest and non integers", { + + ## non integers should give error + expect_error( + temp5 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + boot = FALSE, test = "LRT") + ) + + # Warning with Wald, *non-robust*, boot = F + expect_warning( + temp6 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + boot = FALSE, + test = "Wald", + robust = FALSE) + ) + + # Everything goes through with Wald, robust, boot = F + expect_silent( + temp7 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + boot = FALSE, + test = "Wald", + robust = TRUE) + ) + + # Expect message with noninteger data with boot = T + # because bootstrap is only parametric, and doesn't make sense with EE + expect_message( + # using invisible and capture.output to avoid printed output giving error from a failed + # run of testing testing function + invisible(utils::capture.output(temp8 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + boot = TRUE, + B = 5, + test = "Wald", + robust = TRUE, + allow_noninteger = TRUE))) + ) + + expect_silent( + temp11 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "Wald", + robust=TRUE) + ) + expect_warning( + temp12 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "Wald", + robust=FALSE) + ) + + # Amy needs to think about whether this is allowed!! + expect_error( + temp13 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "LRT", + boot=FALSE) + ) + + expect_error( + temp14 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants + DayAmdmt, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "LRT", + B=3, + boot=TRUE) + ) + + + expect_silent( + temp15 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, + allow_noninteger = TRUE, + test = "Rao", + robust=TRUE, + B=3, + boot=TRUE) + ) + + expect_silent( + temp16 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil, + allow_noninteger = TRUE, + test = "Rao", + robust=FALSE, + B=3, + boot=TRUE) + ) + + expect_silent( + temp17 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "Rao", + robust=TRUE, + B=3, + boot=FALSE) + ) + + # expect warning when non-integers are present and robust = FALSE from Rao test + expect_warning( + temp18 <- differentialTest(formula = ~ Plants + DayAmdmt, + phi.formula = ~ Plants, + formula_null = ~ 1, + phi.formula_null = ~ 1, + data = subsoil_nonint, + allow_noninteger = TRUE, + test = "Rao", + robust=FALSE, + B=3, + boot=FALSE) + ) + + expect_false(identical(temp17$p, temp18$p)) + }) + + +} else { + expect_error(warn_sumexp()) +} From 6d9765803f91b272e4051d10f84a1da1839d2d11 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Sat, 18 Oct 2025 18:25:10 +0200 Subject: [PATCH 2/7] Add clean_taxa_names for SummarizedExperiments --- R/clean_taxa_names.R | 27 ++++++++++++++++++++++++++ R/contrastsTest.R | 2 +- tests/testthat/test_testsRequiringSE.R | 21 ++++++++++++++++++++ 3 files changed, 49 insertions(+), 1 deletion(-) diff --git a/R/clean_taxa_names.R b/R/clean_taxa_names.R index 99c7b47..5bf7fb4 100644 --- a/R/clean_taxa_names.R +++ b/R/clean_taxa_names.R @@ -23,3 +23,30 @@ clean_taxa_names <- function(x, name = "OTU") { warn_phyloseq() } } + +#' Rename taxa +#' +#' Renames taxa to have short human-readable names +#' +#' @param x Object of class \code{SummarizedExperiment} +#' @param name Character, defaults to \code{"OTU"}. Optional. String to use in every taxa name. +#' +#' @details The original taxa names are saved as the \code{original_names} attribute. See the example for an example of how to access the original names. +#' +#' @return Object of class \code{SummarizedExperiment}, with taxa renamed (defaults to OTU1, OTU2, ...), with the original taxa names saved as an attribute. +#' +#' @export +clean_taxa_names_se <- function(x, name = "OTU") { + if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + if (inherits(x, "SummarizedExperiment")) { + attr(x, "original_names") <- row.names(x) + row.names(x) <- paste0(name, seq_len(nrow(x))) + return(x) + } else { + stop("clean_taxa_names_se is intended for SummarizedExperiment objects!") + } + } else { + warn_sumexp() + } + +} diff --git a/R/contrastsTest.R b/R/contrastsTest.R index 13440e4..2383e5e 100644 --- a/R/contrastsTest.R +++ b/R/contrastsTest.R @@ -126,7 +126,7 @@ contrastsTest <- function(formula, phi.formula, M <- rowSums(data) } else { - stop("Input must be either data frame, matrix, or phyloseq object!") + stop("Input must be either data frame, matrix, phyloseq object, or SummarizedExperiment object!") } # Set up output diff --git a/tests/testthat/test_testsRequiringSE.R b/tests/testthat/test_testsRequiringSE.R index 91a0882..1584f83 100644 --- a/tests/testthat/test_testsRequiringSE.R +++ b/tests/testthat/test_testsRequiringSE.R @@ -39,6 +39,27 @@ if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { expect_is(out_phylo, "bbdml") }) + test_that("clean_taxa_names works", { + soil <- soil_phylo_ps %>% + phyloseq::subset_samples(DayAmdmt %in% c(11,21)) + + soil <- SummarizedExperiment::SummarizedExperiment( + colData = data.frame(phyloseq::sample_data(soil)), + assays = list(counts = data.frame(phyloseq::otu_table(soil, taxa_are_rows = TRUE))), + rowData = data.frame(phyloseq::tax_table(soil)) + ) + + tmp1 <- clean_taxa_names_se(soil) + tmp2 <- clean_taxa_names_se(soil, name = "Seq") + + expect_equal(length(tmp1), 7770) + expect_equal(row.names(tmp1[1, ]), "OTU1") + expect_equal(row.names(tmp1[1234, ]), "OTU1234") + expect_equal(row.names(tmp2[2543, ]), "Seq2543") + expect_error(clean_taxa_names_se(c(1,2,3))) + }) + + # tests for contrastsTest() set.seed(1) limma_install <- try(find.package("limma"), silent = TRUE) From 41a11516c2d4806e352c1b29356f4090fa9fe721 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Sat, 18 Oct 2025 18:47:56 +0200 Subject: [PATCH 3/7] Add helper functions for SummarizedExperiments --- DESCRIPTION | 2 +- NAMESPACE | 4 ++++ R/differentialTest.R | 4 ++-- R/genInits.R | 2 +- R/getRestrictionTerms.R | 4 ++-- R/otu_to_taxonomy.R | 30 ++++++++++++++++++++++++++++++ R/plot_differentialTest.R | 15 ++++----------- man/clean_taxa_names_se.Rd | 22 ++++++++++++++++++++++ man/convert_sumexp.Rd | 19 +++++++++++++++++++ man/otu_to_taxonomy_se.Rd | 21 +++++++++++++++++++++ man/warn_sumexp.Rd | 11 +++++++++++ 11 files changed, 117 insertions(+), 17 deletions(-) create mode 100644 man/clean_taxa_names_se.Rd create mode 100644 man/convert_sumexp.Rd create mode 100644 man/otu_to_taxonomy_se.Rd create mode 100644 man/warn_sumexp.Rd diff --git a/DESCRIPTION b/DESCRIPTION index 132c868..c23bf7c 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -17,7 +17,7 @@ License: GPL (>= 2) Imports: stats, utils, VGAM, numDeriv, ggplot2, trust, dplyr, magrittr, detectseparation, scales, rlang Encoding: UTF-8 LazyData: true -RoxygenNote: 7.3.2 +RoxygenNote: 7.3.3 Suggests: knitr, rmarkdown, testthat, diff --git a/NAMESPACE b/NAMESPACE index f26bada..b655ef0 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -12,8 +12,10 @@ export(HDIbetabinom) export(bbdml) export(checkNested) export(clean_taxa_names) +export(clean_taxa_names_se) export(contrastsTest) export(convert_phylo) +export(convert_sumexp) export(coth) export(differentialTest) export(fishZ) @@ -24,6 +26,7 @@ export(invlogit) export(logit) export(lrtest) export(otu_to_taxonomy) +export(otu_to_taxonomy_se) export(pbLRT) export(pbRao) export(pbWald) @@ -35,6 +38,7 @@ export(score) export(waldchisq) export(waldt) export(warn_phyloseq) +export(warn_sumexp) importFrom(magrittr,"%>%") importFrom(rlang,.data) importFrom(stats,simulate) diff --git a/R/differentialTest.R b/R/differentialTest.R index bef29a5..d95a06f 100644 --- a/R/differentialTest.R +++ b/R/differentialTest.R @@ -98,8 +98,8 @@ differentialTest <- function(formula, phi.formula, } else if (inherits(data, "SummarizedExperiment")) { if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { # Set up response - taxanames <- row.names(SummarizedExperiment::rowData(data)) - sample_data <- data.frame(SummarizedExperiment::colData(data)) + taxanames <- row.names(data) + sample_data <- SummarizedExperiment::colData(data) } else { warn_sumexp() } diff --git a/R/genInits.R b/R/genInits.R index e742c1f..d6fd961 100644 --- a/R/genInits.R +++ b/R/genInits.R @@ -33,7 +33,7 @@ genInits <- function(W, M, nstart = 1, use = TRUE) { - init.glm <- eval(parse(text = paste("quasibinomial(link =", link,")"))) + init.glm <- eval(parse(text = paste("quasibinomial(link =", link, ")"))) tmp <- stats::glm.fit(x = X, y = cbind(W, M - W), family = init.glm) b_init <- stats::coef(tmp) # Just use 0.5 for phi_init diff --git a/R/getRestrictionTerms.R b/R/getRestrictionTerms.R index d0c9af9..2965cae 100644 --- a/R/getRestrictionTerms.R +++ b/R/getRestrictionTerms.R @@ -18,9 +18,9 @@ getRestrictionTerms <- function(mod, mod_null = NULL, restrictions = NULL, restr assigns.mu <- attr(mod$X.mu, "assign") assigns.phi <- attr(mod$X.phi, "assign") - sortInteraction <- function(x){ + sortInteraction <- function(x) { # Sorts interaction effects so they always match - sapply(lapply(strsplit(x,":"), sort), paste, collapse = ":") + sapply(lapply(strsplit(x, ":"), sort), paste, collapse = ":") } if (!is.null(mod_null)) { diff --git a/R/otu_to_taxonomy.R b/R/otu_to_taxonomy.R index 98e5573..e1a5f3c 100644 --- a/R/otu_to_taxonomy.R +++ b/R/otu_to_taxonomy.R @@ -24,3 +24,33 @@ otu_to_taxonomy <- function(OTU, data, level = NULL) { warn_phyloseq() } } + + +#' Transform OTUs to their taxonomic label +#' +#' @param OTU String vector. Names of OTU labels in \code{data} +#' @param data \code{phyloseq} object with a taxonomy table +#' @param level (Optional). Character vector. Desired taxonomic levels for output. +#' +#' @importFrom magrittr %>% +#' +#' @return String vector. Names of taxonomic labels matching labels of \code{OTU}. +#' +#' @export +otu_to_taxonomy_se <- function(OTU, data, level = NULL) { + if (!inherits(data, "SummarizedExperiment")) { + stop("This function currently only works for SummarizedExperiment objects.") + } + + if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { + if (!(nrow(SummarizedExperiment::rowData(data)) == 0)) { + if (is.null(level)) { + return(apply(SummarizedExperiment::rowData(data)[OTU,], 1, function(x) {paste(stats::na.omit(x), collapse = "_")})) + } else { + return(sapply(SummarizedExperiment::rowData(data)[OTU, level], function(x) {paste(stats::na.omit(x), collapse = "_")})) + } + } + } else { + warn_sumexp() + } +} diff --git a/R/plot_differentialTest.R b/R/plot_differentialTest.R index 1553ccc..5d18995 100644 --- a/R/plot_differentialTest.R +++ b/R/plot_differentialTest.R @@ -37,17 +37,10 @@ plot.differentialTest <- function(x, level = NULL, data_only = FALSE, ...) { } if (inherits(x$data, "SummarizedExperiment")) { - if (!(nrow(SummarizedExperiment::rowData(x$data)) == 0)) { - if (is.null(level)) { - signif_taxa = apply(SummarizedExperiment::rowData(x$data)[signif_taxa,], 1, function(x) {paste(stats::na.omit(x), collapse = '_')}) - } else { - signif_taxa = sapply(SummarizedExperiment::rowData(x$data)[signif_taxa, level], function(x) {paste(stats::na.omit(x), collapse = '_')}) - } - - if (length(unique(signif_taxa)) != length(unique(x$significant_taxa))) { - # Make sure if repeated taxa add unique otu identifiers - signif_taxa <- paste0(signif_taxa, " (", x$significant_taxa, ")") - } + signif_taxa <- otu_to_taxonomy_se(signif_taxa, x$data, level = level) + if (length(unique(signif_taxa)) != length(unique(x$significant_taxa))) { + # Make sure if repeated taxa add unique otu identifiers + signif_taxa <- paste0(signif_taxa, " (", x$significant_taxa, ")") } } diff --git a/man/clean_taxa_names_se.Rd b/man/clean_taxa_names_se.Rd new file mode 100644 index 0000000..bfb15a6 --- /dev/null +++ b/man/clean_taxa_names_se.Rd @@ -0,0 +1,22 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/clean_taxa_names.R +\name{clean_taxa_names_se} +\alias{clean_taxa_names_se} +\title{Rename taxa} +\usage{ +clean_taxa_names_se(x, name = "OTU") +} +\arguments{ +\item{x}{Object of class \code{SummarizedExperiment}} + +\item{name}{Character, defaults to \code{"OTU"}. Optional. String to use in every taxa name.} +} +\value{ +Object of class \code{SummarizedExperiment}, with taxa renamed (defaults to OTU1, OTU2, ...), with the original taxa names saved as an attribute. +} +\description{ +Renames taxa to have short human-readable names +} +\details{ +The original taxa names are saved as the \code{original_names} attribute. See the example for an example of how to access the original names. +} diff --git a/man/convert_sumexp.Rd b/man/convert_sumexp.Rd new file mode 100644 index 0000000..14a21e5 --- /dev/null +++ b/man/convert_sumexp.Rd @@ -0,0 +1,19 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/convert_sumexp.R +\name{convert_sumexp} +\alias{convert_sumexp} +\title{Function to subset and convert SummarizedExperiment data} +\usage{ +convert_sumexp(data, select) +} +\arguments{ +\item{data}{a \code{SummarizedExperiment} object} + +\item{select}{Name of OTU or taxa to select, must match taxa name in \code{data}} +} +\value{ +A \code{data.frame} object, with elements \code{W} as the observed counts, \code{M} as the sequencing depth, and the sample data with their original names. +} +\description{ +Function to subset and convert SummarizedExperiment data +} diff --git a/man/otu_to_taxonomy_se.Rd b/man/otu_to_taxonomy_se.Rd new file mode 100644 index 0000000..c1d035f --- /dev/null +++ b/man/otu_to_taxonomy_se.Rd @@ -0,0 +1,21 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/otu_to_taxonomy.R +\name{otu_to_taxonomy_se} +\alias{otu_to_taxonomy_se} +\title{Transform OTUs to their taxonomic label} +\usage{ +otu_to_taxonomy_se(OTU, data, level = NULL) +} +\arguments{ +\item{OTU}{String vector. Names of OTU labels in \code{data}} + +\item{data}{\code{phyloseq} object with a taxonomy table} + +\item{level}{(Optional). Character vector. Desired taxonomic levels for output.} +} +\value{ +String vector. Names of taxonomic labels matching labels of \code{OTU}. +} +\description{ +Transform OTUs to their taxonomic label +} diff --git a/man/warn_sumexp.Rd b/man/warn_sumexp.Rd new file mode 100644 index 0000000..ae84cd4 --- /dev/null +++ b/man/warn_sumexp.Rd @@ -0,0 +1,11 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/warn_phyloseq.R +\name{warn_sumexp} +\alias{warn_sumexp} +\title{Function to throw error if the `SummarizedExperiment` package is called but it is not installed} +\usage{ +warn_sumexp() +} +\description{ +Function to throw error if the `SummarizedExperiment` package is called but it is not installed +} From 0b1cdf7ca293ed795fc175899ac430f133353827 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Sat, 18 Oct 2025 18:57:01 +0200 Subject: [PATCH 4/7] Add SummarizedExperiment information to bbdml and differentialTest --- R/bbdml.R | 2 +- R/differentialTest.R | 2 +- man/bbdml.Rd | 2 +- man/differentialTest.Rd | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/R/bbdml.R b/R/bbdml.R index df1089c..ebdc47a 100644 --- a/R/bbdml.R +++ b/R/bbdml.R @@ -2,7 +2,7 @@ #' #' @param formula an object of class \code{formula}: a symbolic description of the model to be fitted to the abundance #' @param phi.formula an object of class \code{formula} without the response: a symbolic description of the model to be fitted to the dispersion -#' @param data a data frame or \code{phyloseq} object containing the variables in the models +#' @param data a data frame, \code{phyloseq}, or \code{SummarizedExperiment} object containing the variables in the models #' @param link link function for abundance covariates, defaults to \code{"logit"} #' @param phi.link link function for dispersion covariates, defaults to \code{"logit"} #' @param method optimization method, defaults to \code{"trust"}, or see \code{\link[optimx]{optimr}} for other options diff --git a/R/differentialTest.R b/R/differentialTest.R index d95a06f..361a5ca 100644 --- a/R/differentialTest.R +++ b/R/differentialTest.R @@ -4,7 +4,7 @@ #' @param phi.formula an object of class \code{formula} without the response: a symbolic description of the model to be fitted to the dispersion #' @param formula_null Formula for mean under null, without response #' @param phi.formula_null Formula for overdispersion under null, without response -#' @param data a data frame containing the OTU table, or \code{phyloseq} object containing the variables in the models +#' @param data a data frame containing the OTU table, \code{phyloseq}, or \code{SummarizedExperiment} object containing the variables in the models #' @param link link function for abundance covariates, defaults to \code{"logit"} #' @param phi.link link function for dispersion covariates, defaults to \code{"logit"} #' @param test Character. Hypothesis testing procedure to use. One of \code{"Wald"}, \code{"LRT"} (likelihood ratio test), or \code{"Rao"}. diff --git a/man/bbdml.Rd b/man/bbdml.Rd index d77e898..5f78b73 100644 --- a/man/bbdml.Rd +++ b/man/bbdml.Rd @@ -25,7 +25,7 @@ bbdml( \item{phi.formula}{an object of class \code{formula} without the response: a symbolic description of the model to be fitted to the dispersion} -\item{data}{a data frame or \code{phyloseq} object containing the variables in the models} +\item{data}{a data frame, \code{phyloseq}, or \code{SummarizedExperiment} object containing the variables in the models} \item{link}{link function for abundance covariates, defaults to \code{"logit"}} diff --git a/man/differentialTest.Rd b/man/differentialTest.Rd index 29a2483..ab1a810 100644 --- a/man/differentialTest.Rd +++ b/man/differentialTest.Rd @@ -38,7 +38,7 @@ differentialTest( \item{phi.formula_null}{Formula for overdispersion under null, without response} -\item{data}{a data frame containing the OTU table, or \code{phyloseq} object containing the variables in the models} +\item{data}{a data frame containing the OTU table, \code{phyloseq}, or \code{SummarizedExperiment} object containing the variables in the models} \item{link}{link function for abundance covariates, defaults to \code{"logit"}} From 334693977fcfe7d2c23104b05677f638a05a2165 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Mon, 20 Oct 2025 20:11:22 +0200 Subject: [PATCH 5/7] Change significant taxa printing only if taxonomic information exists --- R/contrastsTest.R | 2 +- R/otu_to_taxonomy.R | 16 +++++++--------- R/plot_differentialTest.R | 10 ++++++---- 3 files changed, 14 insertions(+), 14 deletions(-) diff --git a/R/contrastsTest.R b/R/contrastsTest.R index 2383e5e..a01dc9b 100644 --- a/R/contrastsTest.R +++ b/R/contrastsTest.R @@ -98,7 +98,7 @@ contrastsTest <- function(formula, phi.formula, } else if (inherits(data, "SummarizedExperiment")) { if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { # Set up response - taxanames <- row.names(SummarizedExperiment::rowData(data)) + taxanames <- row.names(data) sample_data <- SummarizedExperiment::colData(data) } else { warn_sumexp() diff --git a/R/otu_to_taxonomy.R b/R/otu_to_taxonomy.R index e1a5f3c..91e14f4 100644 --- a/R/otu_to_taxonomy.R +++ b/R/otu_to_taxonomy.R @@ -42,15 +42,13 @@ otu_to_taxonomy_se <- function(OTU, data, level = NULL) { stop("This function currently only works for SummarizedExperiment objects.") } - if (requireNamespace("SummarizedExperiment", quietly = TRUE)) { - if (!(nrow(SummarizedExperiment::rowData(data)) == 0)) { - if (is.null(level)) { - return(apply(SummarizedExperiment::rowData(data)[OTU,], 1, function(x) {paste(stats::na.omit(x), collapse = "_")})) - } else { - return(sapply(SummarizedExperiment::rowData(data)[OTU, level], function(x) {paste(stats::na.omit(x), collapse = "_")})) - } - } - } else { + if (!requireNamespace("SummarizedExperiment", quietly = TRUE)) { warn_sumexp() } + + if (is.null(level)) { + apply(SummarizedExperiment::rowData(data)[OTU,], 1, function(x) {paste(stats::na.omit(x), collapse = "_")}) + } else { + sapply(SummarizedExperiment::rowData(data)[OTU, level], function(x) {paste(stats::na.omit(x), collapse = "_")}) + } } diff --git a/R/plot_differentialTest.R b/R/plot_differentialTest.R index 5d18995..f595b86 100644 --- a/R/plot_differentialTest.R +++ b/R/plot_differentialTest.R @@ -37,10 +37,12 @@ plot.differentialTest <- function(x, level = NULL, data_only = FALSE, ...) { } if (inherits(x$data, "SummarizedExperiment")) { - signif_taxa <- otu_to_taxonomy_se(signif_taxa, x$data, level = level) - if (length(unique(signif_taxa)) != length(unique(x$significant_taxa))) { - # Make sure if repeated taxa add unique otu identifiers - signif_taxa <- paste0(signif_taxa, " (", x$significant_taxa, ")") + if (!(nrow(SummarizedExperiment::rowData(x$data)) == 0) && "kingdom" %in% tolower(colnames(SummarizedExperiment::rowData(x$data)))) { + signif_taxa <- otu_to_taxonomy_se(signif_taxa, x$data, level = level) + if (length(unique(signif_taxa)) != length(unique(x$significant_taxa))) { + # Make sure if repeated taxa add unique otu identifiers + signif_taxa <- paste0(signif_taxa, " (", x$significant_taxa, ")") + } } } From 00c754a578632f3850678566c0312b210f846159 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Fri, 31 Oct 2025 21:22:13 +0200 Subject: [PATCH 6/7] Add SummarizedExperiment to github workflows --- .github/workflows/R-CMD-check.yaml | 2 +- .github/workflows/test-coverage.yaml | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index 450b1f4..8f6c68a 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -64,7 +64,7 @@ jobs: run: | remotes::install_deps(dependencies = TRUE) remotes::install_cran("rcmdcheck") - if (!requireNamespace("BiocManager", quietly = TRUE)){install.packages("BiocManager")}; BiocManager::install(c("phyloseq", "limma"), ask = FALSE) + if (!requireNamespace("BiocManager", quietly = TRUE)){install.packages("BiocManager")}; BiocManager::install(c("phyloseq", "limma", "SummarizedExperiment"), ask = FALSE) shell: Rscript {0} - name: Check diff --git a/.github/workflows/test-coverage.yaml b/.github/workflows/test-coverage.yaml index c62d30c..0fb201b 100644 --- a/.github/workflows/test-coverage.yaml +++ b/.github/workflows/test-coverage.yaml @@ -37,7 +37,7 @@ jobs: - name: Install dependencies run: | install.packages(c("remotes")) - if (!requireNamespace("BiocManager", quietly = TRUE)){install.packages("BiocManager")}; BiocManager::install(c("phyloseq", "limma"), ask = FALSE) + if (!requireNamespace("BiocManager", quietly = TRUE)){install.packages("BiocManager")}; BiocManager::install(c("phyloseq", "limma", "SummarizedExperiment"), ask = FALSE) remotes::install_deps(dependencies = TRUE) remotes::install_cran("covr") shell: Rscript {0} From 6572a52880aed485f81fe3f05f84a6cdc3192aa5 Mon Sep 17 00:00:00 2001 From: Luke Hannan Date: Fri, 31 Oct 2025 21:42:32 +0200 Subject: [PATCH 7/7] Add SummarizedExperiment information to vignette --- vignettes/corncob-intro.Rmd | 3 +++ 1 file changed, 3 insertions(+) diff --git a/vignettes/corncob-intro.Rmd b/vignettes/corncob-intro.Rmd index 8252479..07175ff 100644 --- a/vignettes/corncob-intro.Rmd +++ b/vignettes/corncob-intro.Rmd @@ -41,6 +41,9 @@ knitr::opts_chunk$set(fig.width=8, fig.height=4) Here, we introduce `corncob`, an individual taxon regression model that uses abundance tables and sample data. `corncob` is able to model differential abundance and differential variability, and addresses each of the challenges presented above. Note that in order to follow along with this tutorial (but not to use `corncob`!) you will need to have `phyloseq` installed. We will check if you have `phyloseq` installed, and if you do not then you can read the following code but it will not be run. See the vignette `corncob-intro-no-phyloseq.Rmd` for a version of this vignette without a dependence on `phyloseq`. + +If your data is in a `TreeSummarizedExperiment` object, or any object that inherits from `SummarizedExperiment` `corncob` will still work. `corncob` will assume that the first assay in the assay list of your `SummarizedExperiment` object is the abundance table of interest, and that the `colData` slot contains information about the samples. + ```{r, results = 'hide'} phy <- requireNamespace("phyloseq", quietly = TRUE) == TRUE ```