---
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()
```
