--- title: "DAssemble: Ensemble Differential Analysis" author: - name: Himel Mallick - name: Ziyu Liu - name: Nalin Arora package: DAssemble output: BiocStyle::html_document: toc: true toc_float: true vignette: > %\VignetteIndexEntry{DAssemble: Ensemble Differential Analysis} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include=FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", warning = FALSE ) ``` # Introduction DAssemble is an R package for ensemble-based differential association analysis of high-throughput omics data. It provides a flexible framework for combining a user-specified **core method** with one or more complementary **enhancer** methods, allowing users to integrate statistical evidence across different modeling assumptions rather than relying on a single differential analysis approach. DAssemble aggregates model-specific p-values using the **Cauchy Combination Tests (CCT)** and reports both ensemble-level and method-specific results. The package was developed to support robust and interpretable differential analysis across multiple omics domains, including bulk RNA-seq, single-cell RNA-seq, and microbiome data. We submitted DAssemble to Bioconductor to make the method easier to install, document, test, and integrate with established Bioconductor workflows and data structures, thereby supporting reproducible and extensible analysis of high-throughput biological data. # Overview `DAssemble` implements an ensemble framework for differential-abundance and differential-expression analysis. Users can choose a single core model and optionally combine it with enhancer tests using the Cauchy Combination Test. The package accepts a `MultiAssayExperiment` object directly. When multiple experiments are present, use `assay_name` to select the experiment to analyze. For lower-level workflows, `features` can also be a data frame with samples in rows and features in columns, paired with a `metadata` data frame with matching sample row names. The exposure variable supplied through `expVar` must be binary. # Installation Install `DAssemble` from Bioconductor with `BiocManager`: ```{r installation, eval=FALSE} if (!requireNamespace("BiocManager", quietly = TRUE)) { install.packages("BiocManager") } BiocManager::install("DAssemble") ``` # Bulk RNA-seq Example The `airway` dataset from Bioconductor contains RNA-seq counts from airway smooth muscle cells treated with dexamethasone. The example below runs `DAssemble` with `DESeq2` as the core method and two enhancers. ```{r airway-example} data("airway", package = "airway") counts <- SummarizedExperiment::assay(airway, "counts") keep_genes <- order(rowSums(counts), decreasing = TRUE)[seq_len(500L)] counts <- counts[keep_genes, , drop = FALSE] metadata <- as.data.frame(SummarizedExperiment::colData(airway)) metadata <- metadata[colnames(counts), , drop = FALSE] se <- SummarizedExperiment::SummarizedExperiment( assays = list(counts = counts) ) mae <- MultiAssayExperiment::MultiAssayExperiment( experiments = list(rnaseq = se), colData = S4Vectors::DataFrame(metadata) ) res <- DAssemble::DAssemble( features = mae, assay_name = "rnaseq", core_method = "DESeq2", enhancers = c("WLX", "LR"), expVar = "dex", p_adj = "BH", enhancer_norm = "tmm", return_components = TRUE, return_subensembles = TRUE ) head(res$res) names(res$components) ``` This analysis combines DESeq2 p-values with Wilcoxon and logistic-regression enhancer p-values. The combined results are returned in `res$res`, and the per-method outputs are available in `res$components`. # Microbiome Example The `GlobalPatterns` dataset from `phyloseq` contains 16S rRNA profiles from environmental and host-associated microbiome samples. The example below compares fecal and soil samples using enhancer-only DAssemble with CLR normalization. ```{r microbiome-example} data("GlobalPatterns", package = "phyloseq") gp <- GlobalPatterns otu <- as(phyloseq::otu_table(gp), "matrix") metadata <- as.data.frame(phyloseq::sample_data(gp)) keep <- metadata$SampleType %in% c("Feces", "Soil") otu <- otu[, keep, drop = FALSE] keep_taxa <- order(rowSums(otu), decreasing = TRUE)[seq_len(250L)] features <- as.data.frame(t(otu[keep_taxa, , drop = FALSE])) metadata <- metadata[keep, , drop = FALSE] metadata$group <- droplevels(factor(metadata$SampleType)) stopifnot(nlevels(metadata$group) == 2L) se <- SummarizedExperiment::SummarizedExperiment( assays = list(counts = t(as.matrix(features))) ) mae <- MultiAssayExperiment::MultiAssayExperiment( experiments = list(microbiome = se), colData = S4Vectors::DataFrame(metadata) ) res <- DAssemble::DAssemble( features = mae, assay_name = "microbiome", core_method = NULL, enhancers = c("WLX", "LR", "KS"), expVar = "group", enhancer_norm = "clr", return_components = TRUE, return_subensembles = TRUE ) head(res$res) names(res$components) head(res$ensembles) ``` Because the core method is `NULL`, this example combines only enhancer tests. The `$ensembles` element reports the CCT results for the enhancer combinations when `return_subensembles = TRUE`. # Single-Cell Example The single-cell example uses a small `SingleCellExperiment` object with two cell groups. The object is wrapped in a `MultiAssayExperiment` before being passed to `DAssemble`, matching the same input architecture used above. ```{r single-cell-example} set.seed(1) n_genes <- 120L n_cells <- 20L cell_type <- factor(rep(c("control", "treated"), each = n_cells / 2L)) counts <- matrix( stats::rnbinom(n_genes * n_cells, mu = 20, size = 5), nrow = n_genes, dimnames = list( paste0("gene", seq_len(n_genes)), paste0("cell", seq_len(n_cells)) ) ) counts[seq_len(12L), cell_type == "treated"] <- counts[seq_len(12L), cell_type == "treated"] + 15L metadata <- data.frame( CellType = cell_type, row.names = colnames(counts) ) sce <- SingleCellExperiment::SingleCellExperiment( assays = list(counts = counts), colData = S4Vectors::DataFrame(metadata) ) mae <- MultiAssayExperiment::MultiAssayExperiment( experiments = list(single_cell = sce), colData = S4Vectors::DataFrame(metadata) ) res <- DAssemble::DAssemble( features = mae, assay_name = "single_cell", core_method = "edgeR", enhancers = c("WLX", "LR", "KS"), enhancer_norm = "TSS", expVar = "CellType", p_adj = "BH", return_components = TRUE, return_subensembles = TRUE ) head(res$res) names(res$components) ``` # Session information ```{r session-info} sessionInfo() ```