Package {ALDEx3}

ALDEx3 Linear Models

Description

ALDEx3 Linear Models

Usage

aldex(
  Y,
  X,
  data = NULL,
  method = "lm",
  nsample = 2000,
  scale = NULL,
  streamsize = 8000,
  n.cores = detectCores() - 1,
  return.pars = c("X", "estimate", "std.error", "p.val", "p.val.adj", "logComp",
    "logScale"),
  p.adjust.method = "BH",
  test = "t.HC3",
  onesided = FALSE,
  ...
)

Arguments

Details

Sequence counts determine the relative composition of each sample, but not its total abundance or scale. The scale argument tells ALDEx3 what is known about that missing quantity and how uncertain that knowledge is.

Choose a scale model based on the information available:

• Use sample.sm for external measurements such as flow cytometry, qPCR, or spike-ins for individual samples.

• Use coefficient.sm when prior knowledge concerns a group difference, treatment effect, time trend, or another model term.

• Use clr.sm when CLR normalization is a reasonable center but its implied scale differences are uncertain.

• Use tss.sm when the starting assumption is no systematic scale difference between groups.

• Write a custom scale-model function when the built-in models do not match the information available for the study. The "Writing a custom scale model" example below describes the required function interface and N x nsample output.

Value

a list with elements controled by parameter resturn.pars. Options include: - X: P x N covariate matrix - estimate: (P x D x nsample) array of linear model estimates - std.error: (P x D x nsample) array of standard error for the estimates - p.val: (P x D) matrix, unadjusted p-value for two-sided t-test - p.val.adj: (P x D) matrix, p-value for two-sided t-test adjusted according to p.adj.method - logScale: (N x S) matrix, samples of the log scale from the scale model - logComp: (D x N x S) array, samples of the log composition from the multinomial-Dirichlet - streaming: boolean, detnote if streaming was used. - random.effects (Pr x N x S): if using mixed effects models, return all Pr random-effects covariance terms and the residual variance. For method = "blmm", correlated random-effect covariance terms are returned when present, matching the lme4 naming convention as closely as possible. Note, logScale and logComp are not returned if streaming is active.

Author(s)

Justin Silverman, Kyle McGovern

Examples

Y <- matrix(1:110, 10, 11)
condition <- c(rep(0, 5), rep(1, 6))
data <- data.frame(condition=condition)
## Use CLR as the center and allow coefficient-scale uncertainty of SD 0.5.
res <- aldex(Y, ~condition, data, nsample=2000, scale=clr.sm, gamma=0.5)

## Writing a custom scale model
## A custom model must return one log2-scale draw per sample and Monte
## Carlo replicate. Here, external_scale supplies a known value per sample.
known_scale <- function(X, logComp, external_scale) {
  nsample <- dim(logComp)[3]
  matrix(rep(external_scale, nsample), nrow=ncol(X), ncol=nsample)
}
res <- aldex(Y, ~condition, data, nsample=2000, scale=known_scale,
             external_scale=rep(0, ncol(Y)))

ALDEx3 Effect Diagnostics

Description

ALDEx2-inspired effect diagnostics for an ALDEx3 binary contrast

Usage

aldex.effect(object, contrast)

Arguments

Details

aldex.effect is a helper for effect-size diagnostics that Greg Gloor and ALDEx2 users commonly inspect when evaluating pairwise effects. It is not a replacement for summary.aldex, which reports model estimates, standard errors and adjusted p-values. This function summarizes a single binary ALDEx3 model contrast using Cohen's-d-based Monte Carlo diagnostics.

For a selected contrast row k, feature d, Monte Carlo draw s, and binary group indicators x = object$X[k, ], let B[d, s] be the fitted ALDEx3 coefficient stored in object$estimate[k, d, s]. Let logW[d, i, s] = object$logComp[d, i, s] + object$logScale[i, s] be the reconstructed log abundance for sample i. aldex.effect requires logComp and logScale, so it is unavailable when aldex streams large jobs and omits those arrays. For each feature and Monte Carlo draw, the pooled within-group variance is ((n0 - 1) * var(logW[d, x == 0, s]) + (n1 - 1) * var(logW[d, x == 1, s])) / (n0 + n1 - 2), where n0 and n1 are the group sizes.

The returned columns are calculated as follows:

• estimate: the mean of B[d, s] over Monte Carlo draws.

• pooled.SD: the mean of the pooled within-group standard deviation, ⁠sqrt(pooled variance)⁠, over Monte Carlo draws.

• cohens.d: the mean of B[d, s] / pooled.SD[d, s] over Monte Carlo draws. This is Cohen's-d-based and is not the same as the historical robust ALDEx2 effect statistic based on between/within dispersion.

• overlap: an ALDEx2-style directional uncertainty diagnostic. It counts the Monte Carlo Cohen's d draws below and above zero, smooths those two counts with a 0.5 pseudocount, converts them with the ALDEx/Aitchison mean, and returns the smaller smoothed sign probability. Values near 0 indicate that nearly all Monte Carlo effect draws have the same sign. Values near 0.5 indicate uncertainty about direction. overlap is not a p-value, a confidence interval probability, or a literal area of overlap between two density curves.

Value

A data.frame with columns parameter, entity, estimate, pooled.SD, cohens.d, and overlap.

Author(s)

Greg Gloor, Justin Silverman

Simulation function soley for testing and exploring ALDEx3, Truth is in CLR Coordinates.

Description

Not designed to create realistic data. Does not add any noise to linear regression! True W is in CLR Corrdinates

Usage

aldex.lm.sim.clr(D = 10, N = 11, P = 2, depth = 10000)

Arguments

Value

a list with elements Y, X, W, and Lambda

Author(s)

Justin Silverman

Simulation for testing mixed effects models.

Description

Includes random intercpt, option to include random slope, second random intercept, and time-correlation.

Usage

aldex.mem.sim(
  D,
  days,
  subjects,
  depth = 10000,
  location = 1,
  random_slope = FALSE,
  corr = 0,
  rho_ar1 = 0,
  sd_resid = 0.1
)

Arguments

Author(s)

Kyle McGovern

Plot Method for ALDEx3 objects with pairwise comparisons

Description

Simple plots for an ALDEx3 object

Usage

aldex.plot(
  object,
  plot = c("volcano", "effect", "MA", "water"),
  contrast = NULL,
  threshold = 0.05,
  min.diff = 0.5,
  cohen = 0.5,
  sig.col = rgb(1, 0, 0, 0.5),
  water.show = 5,
  water.col = c("red", "blue"),
  water.names = TRUE,
  ...
)

Arguments

Details

Provides volcano, effect, MA and waterfall plots from an ALDEx3 object.

This method plots combinations of adjusted p-values from object$p.val.adj, posterior estimates, standard errors and log abundance values averaged across Monte Carlo samples. The result is returned as a single plot of the desired type. Effect plots use the ALDEx2-inspired diagnostics returned by aldex.effect. The contrast must identify exactly one binary model coefficient. Effect diagnostics and MA plots require logComp and logScale, which are not returned when aldex streams large jobs.

Value

the desired plot

Author(s)

Greg Gloor

Calculate p-values adjusting for changes in sign as described by Nixon et al. (2024) in Beyond Normalization: Incorporating Scale Uncertainty in Microbiome and Gene Expression Analysis (internal only)

Description

Calculate p-values adjusting for changes in sign as described by Nixon et al. (2024) in Beyond Normalization: Incorporating Scale Uncertainty in Microbiome and Gene Expression Analysis (internal only)

Usage

aldex.pvals(p.lower, p.upper, p.adjust.method, onesided = FALSE)

Arguments

Value

A list with P x D matrices with the non-adjusted and adjusted p-values.

Author(s)

Justin Silverman, Kyle McGovern

References

Nixon G, Gloor GB, Silverman JD (2025). "Incorporating scale uncertainty in microbiome and gene expression analysis as an extension of normalization". Genome Biology. doi:10.1186/s13059-025-03609-3

Fast approximate mixed-effects engine for ALDEx3.

Description

Fast approximate mixed-effects engine for ALDEx3.

Usage

blmm(logW, formula, data, n.cores = 1L)

Arguments

Details

Internal engine called by aldex when method = "blmm".

Uses an lme4-style profiled Gaussian LMM formulation. One batched anchor REML problem is fit per feature across all S MC draws; draw-specific local covariance updates are obtained via a single Newton step; exact conditional GLS fixed-effect solves are then run per draw. The approximation is only in the variance-component step — fixed effects and standard errors are exact conditional on the updated covariance parameters.

Features for which the approximate path fails are silently re-fit with the exact lme4 engine and a consolidated warning is issued naming the affected features and their error messages.

Value

named list of arrays (P x D x S): estimate, std.error, df, p.lower, p.upper, and (Pr x D x S) random.eff — same layout as sr.mem.

CLR-based scale model

Description

Use this model when centered log-ratio (CLR) normalization is a reasonable starting point, but you do not want to assume that it recovers the sample scales exactly. For example, in a case-control study, CLR may provide a plausible center while gamma allows the true difference in total abundance between cases and controls to depart from that center.

Usage

clr.sm(X, logComp, gamma = 0.5)

Arguments

Details

Set gamma = 0 to use the CLR-implied scales without additional uncertainty. A value such as gamma = 0.5 allows each scale-model coefficient to vary with SD 0.5 on the log2 scale. Larger values express less confidence in the CLR assumption.

The uncertainty applies to model coefficients, not independently to every sample. Samples with the same covariate values therefore receive the same random shift in a Monte Carlo draw. This differs from sample.sm, which can draw a separate scale for each sample.

Value

A numeric matrix of dimension N x nsample giving Monte Carlo samples of the log-scale for each sample (rows) and each Monte Carlo draw (columns).

Author(s)

Justin Silverman

References

Nixon G, Gloor GB, Silverman JD (2025). "Incorporating scale uncertainty in microbiome and gene expression analysis as an extension of normalization". Genome Biology. doi:10.1186/s13059-025-03609-3

Examples

# Allow moderate uncertainty around CLR-implied scale differences.
Y <- matrix(seq_len(60), nrow = 10)
condition <- factor(c("control", "control", "control",
                      "case", "case", "case"))
fit <- aldex(Y, ~ condition, data.frame(condition), nsample = 20,
             scale = clr.sm, gamma = 0.5)

Scale model using prior information about groups or treatments

Description

Use this model when prior information describes a treatment, group, time, or other model coefficient rather than an individual sample. Samples with the same covariate values share the same scale shift in each Monte Carlo draw.

Usage

coefficient.sm(
  X,
  logComp,
  c.mu = NULL,
  c.cor = NULL,
  c.sd = NULL,
  c.cov = NULL
)

Arguments

Details

Suppose an effect plot from a control-versus-treatment experiment is centered near -8, even though prior knowledge suggests that most features should not change. One way to represent that knowledge is to shift the treatment samples upward by 8 log2 units relative to the controls. With control as the reference level, use c.mu = c(0, 8): the first value leaves the common intercept unchanged, and the second adds 8 to treatment. Setting c.sd = c(0, 0.5) expresses an SD of 0.5 around that treatment shift. If the plot were centered near +8, use -8 instead. Adding 8 to both groups would not recenter the plot because it would leave their difference unchanged.

More generally, c.mu gives the expected log2 scale change for each model coefficient and c.sd gives its SD. Use c.cov when the coefficient uncertainties are correlated. For each Monte Carlo draw, ALDEx3 samples a coefficient vector and applies it to the model matrix X; mathematically, b^{(m)} \sim N(c.mu, Sigma) and the sample scales are b^{(m)T} X.

Exactly one of c.sd or c.cov must be provided. The deprecated c.cor argument may be used in place of c.cov during the compatibility period.

Value

A numeric matrix of dimension N x nsample giving Monte Carlo draws of the log2 scale for each sample (rows) across nsample draws (columns).

Author(s)

Kyle McGovern

Examples

# Three control samples followed by three treatment samples.
condition <- factor(
  c("control", "control", "control", "treatment", "treatment", "treatment"),
  levels = c("control", "treatment")
)
X <- t(model.matrix(~ condition))
logComp <- array(0, c(2, 6, 10))
logScale <- coefficient.sm(
  # Shift treatment upward by 8, with SD 0.5 around that shift.
  X, logComp, c.mu = c(0, 8), c.sd = c(0, 0.5)
)

Cohen's D

Description

Function to compute cohensd on the results provided by the aldex function

Usage

cohensd(m, var)

Arguments

Details

WARNING: this function is experimental and requires users read the documentation fully.

Value

A (D x nsample)-matrix of Cohen's d statistics for the variable of interest.

Author(s)

Justin Silverman

Combine output from aldex streams (internal only)

Description

Takes a list l, where every element is itself a list of 3D arrays. Each array should have matching first and second dimentions. This function combines those arrays along the third dimension.

Usage

combine.streams(l)

Arguments

Value

a list of 3D arrays

Author(s)

Justin Silverman

Freaking Fask Linear Models

Description

Tailored for ALDEx3 where covariates are shared between massive numbers of linear regressions where only Y is changing.

Usage

fflm(Y, X, test = "t.HC3")

Arguments

Value

A list of (P x D x S)-arrays with the OLS point estimates, the standard errors, and the two-sided p-values for each coefficient (P), of each model fit to each taxa (D) and each posterior sample (S)

Author(s)

Justin Silverman

Gut Crohn's microbiome dataset (list)

Description

This dataset contains a matrix and a data frame: genus-level microbiome profiles and corresponding sample metadata from a Crohn's disease case–control cohort. The dataset is used in examples and vignettes throughout the package.

Usage

gut_crohns_data

Format

A named list of one matrix and a dataframe:

counts

A matrix with read counts with 195 rows and 95 columns)

metadata

A data frame with 95 rows and 7 columns containing subject-level covariates.

Details

The counts matrix has one row per sample and one column per genus. The metadata data frame has one row per sample with metadata, critically Health.status either CD or Control, and Average cell count per gram frozen feces.

Source

Vandeputte D, Falony G, Vieira-Silva S, Wang J, Sailer M, Theis S, Raes J (2017). "Quantitative microbiome profiling links gut community variation to microbial load". Nature, 551, 507–511. doi:10.1038/nature24460

Closure operation

Description

Closure operation

Usage

miniclo(X)

Arguments

Value

D x N matrix with columns summing to 1

Author(s)

Justin Silverman

Oral microbiome perturbation dataset (list)

Description

This dataset contains a matrix and a data frame: genus-level microbiome profiles and corresponding sample metadata from an oral microbiome perturbation study. 28 participants' oral microbiomes were measured before, 15 minutes after, and 2 hours after perturbation with either a water control, antiseptic mouthwash, alchohol-free mouthwash, or soda. The dataset is used in examples and vignettes throughout the package.

Usage

oral_mouthwash_data

Format

A named list of one matrix and a dataframe:

counts

A matrix with read counts with 116 rows and 81 columns)

metadata

A data frame with 81 rows and 38 columns containing sample-level covariates.

Details

The counts matrix has one row per sample and one column per genus. The metadata data frame has one row per sample with metadata, critically the participant ID, flow cytometry average (and replicate) cells, time_c (time points), and treat (the perturbations )

Source

Marotz C, Morton JT, Navarro P, Coker J, Belda-Ferre P, Knight R (2021). "Quantifying live microbial load in human saliva samples over time reveals stable composition and dynamic load". mSystems, 6(3), e01182-21. doi:10.1128/mSystems.01182-20

Function for sampling Dirichlet random variables (base-2 normalized via log-sum-exp for stability)

Description

Function for sampling Dirichlet random variables (base-2 normalized via log-sum-exp for stability)

Usage

rLogDirichlet(n, alpha)

Arguments

Value

a D x n matrix of samples

Author(s)

Justin Silverman

Test object contains elements

Description

Test object contains elements

Usage

req(obj, names)

Arguments

Author(s)

Justin Silverman

Per-sample scale model using external information

Description

Use this model when each sample has its own external scale measurement, such as a flow-cytometry cell count, qPCR concentration, or spike-in estimate. Put the log2 measurement for each sample in s.mu. Use s.sd to describe the measurement uncertainty for each sample.

Usage

sample.sm(
  X,
  logComp,
  s.mu = NULL,
  s.var = NULL,
  s.cor = NULL,
  s.sd = NULL,
  s.cov = NULL
)

Arguments

Details

For example, if four samples have estimated total cell counts of 1e8, 2e8, 8e7, and 1.5e8, set s.mu = log2(c(1e8, 2e8, 8e7, 1.5e8)). If each measurement has an SD of 0.5 on the log2 scale, set s.sd = rep(0.5, 4).

Use s.sd when measurement errors are independent. Use s.cov instead when errors are correlated, for example because several samples share a calibration standard. The diagonal of s.cov contains variances, and its off-diagonal entries contain covariances. Thus independent SDs of 0.5 can also be written as s.cov = diag(0.5^2, 4).

Supply exactly one of s.sd or s.cov. With s.sd, each sample receives its own independent scale draw. Samples in the same group do not automatically share uncertainty as they do in coefficient-based models.

Value

A numeric matrix of dimension N x nsample giving Monte Carlo draws of the log2 scale for each sample (rows) across nsample draws (columns).

Author(s)

Kyle McGovern

Examples

# External total-cell estimates for four samples.
cell_count <- c(1e8, 2e8, 8e7, 1.5e8)
X <- matrix(1, 1, 4)
logComp <- array(0, c(2, 4, 10))
logScale <- sample.sm(
  X, logComp, s.mu = log2(cell_count), s.sd = rep(0.5, 4)
)

Implementation of SR-MEM: scale-reliant mixed effects models.

Description

Implementation of SR-MEM: scale-reliant mixed effects models.

Usage

sr.mem(logW, formula, data, n.cores, method, mem.args)

Arguments

Value

A list of (P x D x S)-arrays with the fixed effect point estimates, the standard errors, degrees of freedom and the lower and upper p-values for each coefficient (P), of each model fit to each taxa (D) and each posterior sample (S) and a (Pr x D x S)-array with (Pr) random effects.

Author(s)

Kyle McGovern

Summary Method for ALDEx3 Objects

Description

Summarize an ALDEx3 result object

Usage

## S3 method for class 'aldex'
summary(object, ignore.intercept = TRUE, ...)

Arguments

Details

Provides a summary of the adjusted p-values, estimates, and standard errors from an ALDEx3 result object.

This method extracts adjusted p-values from object$p.val.adj, along with posterior estimates and standard errors averaged across Monte Carlo samples. The result is returned as a long-format data.frame suitable for downstream analysis or visualization.

Value

A data.frame with columns parameter, entity, p.val.adjusted, estimate, and std.error.

Author(s)

Justin Silverman

TSS-centered scale model

Description

Total-sum scaling (TSS) starts from the assumption that total scale does not change systematically between experimental groups. Use this model when that is a reasonable starting point but you want to express uncertainty about it. For example, if a treatment is expected to redistribute features without changing total biomass, TSS centers the treatment-versus-control scale difference at zero.

Usage

tss.sm(X, logComp, gamma = 0.5)

Arguments

Details

Set gamma = 0 to enforce no scale difference. Setting gamma = 0.5 keeps the expected difference at zero but allows each scale-model coefficient to vary with SD 0.5 on the log2 scale. Larger values express less confidence that total scale is unchanged.

ALDEx3 applies this uncertainty to model coefficients: for each draw, b^{(m)} \sim N(0, \gamma^2 I), and the sample scales are b^{(m)T} X. Samples with the same covariate values therefore share a random shift. This differs from sample.sm, which can draw independent uncertainty for every sample.

Value

A numeric matrix of dimension N x nsample giving Monte Carlo draws of the log2 scale for each sample (rows) across nsample draws (columns).

Author(s)

Justin Silverman

References

Nixon G, Gloor GB, Silverman JD (2025). "Incorporating scale uncertainty in microbiome and gene expression analysis as an extension of normalization". Genome Biology. doi:10.1186/s13059-025-03609-3

Examples

# Center the group-scale difference at zero, with moderate uncertainty.
Y <- matrix(seq_len(60), nrow = 10)
condition <- factor(c("control", "control", "control",
                      "treatment", "treatment", "treatment"))
fit <- aldex(Y, ~ condition, data.frame(condition), nsample = 20,
             scale = tss.sm, gamma = 0.5)