```{css, echo=FALSE} table td, table th { white-space: nowrap; }
```r knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.width = 7, fig.height = 5 ) library(linf)
The linf package implements L-infinity normalization and Dominant Community
State Types (dCSTs) for compositional data. This vignette demonstrates the
complete workflow in two parts: first a quick-start with toy data to introduce
the core functions, then a real-data analysis of gut microbiome samples from
the American Gut Project showing how to inspect dominant and subdominant
community structure.
L-infinity normalization scales each row of a count matrix by its maximum value, projecting samples onto the L-infinity unit ball. Every normalized row has a maximum of exactly 1.
set.seed(1) S.counts <- matrix( rpois(30, lambda = 5), nrow = 10, ncol = 3, dimnames = list(paste0("s", 1:10), c("Taxon_A", "Taxon_B", "Taxon_C")) ) Z <- normalize.linf(S.counts) apply(Z, 1, max) # all 1
Each sample is assigned to its dominant feature: the column achieving the within-sample maximum. Samples with the same dominant feature form a depth-1 dominance sample set.
cells <- linf.cells(Z) table(cells$label)
linf.csts() groups samples by dominant feature and collapses small groups
(below threshold n0) into a RARE_DOMINANT bucket. This produces
Dominant Community State Types (dCSTs).
A <- matrix(c(5, 1, 0), nrow = 6, ncol = 3, byrow = TRUE) B <- matrix(c(1, 5, 0), nrow = 2, ncol = 3, byrow = TRUE) C <- matrix(c(1, 0, 5), nrow = 2, ncol = 3, byrow = TRUE) S <- rbind(A, B, C) S <- sweep(S, 1, rowSums(S), "/") colnames(S) <- c("Dom1", "Dom2", "Dom3") res <- linf.csts(S, n0 = 5) table(res$cell.label, useNA = "ifany")
Only Dom1 (6 samples) meets the n0 = 5 threshold. Dom2 and Dom3 are
collapsed into RARE_DOMINANT.
Single-species dominance is common in some microbial ecosystems and less common in others. dCSTs provide a deterministic way to summarize the dominant feature and then refine large groups by subdominant features. The gut example below demonstrates those mechanics without treating the selected package subset as an epidemiologic sample.
The agp_gut dataset bundled with this package contains r nrow(agp_gut$counts)
gut microbiome samples from the American Gut Project (PRJEB11419), a citizen-science
16S rRNA V4 survey. The dataset is a dCST-stratified demonstration subset
that includes all samples in four selected uncommon dCSTs and fills
the remaining slots with a seed-42 simple random sample from the eligible
background. Phenotypes are joined only after selection and do not influence
membership. Because inclusion probabilities differ by dCST, it is not a
probability sample of the underlying cohort; phenotype frequencies, effect
sizes, and hypothesis tests must not be interpreted as population results.
data(agp_gut) str(agp_gut, max.level = 1) dim(agp_gut$counts) # samples x taxa
Metadata includes self-reported health conditions parsed from the AGP questionnaire:
disease_cols <- c("IBS", "IBD", "Obesity", "Cardiovascular_disease", "Autoimmune", "Acid_reflux") sapply(disease_cols, function(col) sum(agp_gut$meta[[col]], na.rm = TRUE))
filter.asv() removes low-depth samples and rare taxa in one call. We
require at least 1,000 reads per sample and taxa present in at least 5% of
samples.
filt <- filter.asv(agp_gut$counts, min.lib = 1000, prev.prop = 0.05, min.count = 2) dim(filt$counts) dim(filt$rel)
M <- normalize.linf(filt$counts) # Verify: every row max is 1 stopifnot(all(abs(apply(M, 1, max) - 1) < 1e-10))
csts <- linf.csts(M, n0 = 30) # dCST size distribution sort(table(csts$cell.label), decreasing = TRUE)
The landscape is dominated by Bacteroides and Escherichia-Shigella, with
Prevotella, Faecalibacterium, and several smaller dCSTs forming the tail.
The RARE_DOMINANT bucket collects samples dominated by taxa too infrequent
to form their own dCST at this threshold.
Note on Escherichia-Shigella: this genus is known to be inflated in 16S V4 data due to primer cross-reactivity. Its high prevalence should be interpreted with caution.
refine.linf.csts() subdivides large dCSTs by their subdominant species.
For instance, a Bacteroides-dominated sample might be further classified by
whether its second-most-abundant taxon is Faecalibacterium or
Lachnospiraceae. This captures co-dominance patterns that are especially
important in the diverse gut environment.
csts2 <- refine.linf.csts(M, csts, n0 = 30) # Show depth-2 dCSTs with >= 20 samples tab2 <- sort(table(csts2$cell.label), decreasing = TRUE) tab2[tab2 >= 20]
A useful diagnostic: the relative abundance of the dominant species in each sample. In the gut, most samples have modest dominance (< 50%), consistent with the diverse nature of the ecosystem.
rel <- filt$rel dom_strength <- apply(rel, 1, max) hist(dom_strength, breaks = 50, col = "#2ecc71", border = "white", main = "Dominance Strength in Gut Microbiome", xlab = "Relative abundance of dominant species", ylab = "Number of samples") abline(v = 0.5, col = "red", lty = 2, lwd = 2) legend("topright", "50% dominance", col = "red", lty = 2, lwd = 2, bty = "n")
tab1 <- sort(table(csts$cell.label), decreasing = TRUE) par(mar = c(10, 4, 3, 1)) bp <- barplot(tab1, col = ifelse(names(tab1) == "RARE_DOMINANT", "#e74c3c", "#3498db"), las = 2, cex.names = 0.7, ylab = "Number of samples", main = "Gut dCST Size Distribution (depth 1)") text(bp, tab1 + 5, labels = tab1, cex = 0.7, pos = 3)
The bundle retains selected self-reported phenotype fields so that users can inspect the object structure and practice data alignment. The following table is purely descriptive; dCST-stratified selection precludes population prevalence estimation or treating association calculations as cohort results.
meta <- agp_gut$meta[match(rownames(M), agp_gut$meta$Run), ] stopifnot(all(meta$Run == rownames(M))) conditions <- c("IBS", "IBD", "Obesity", "Cardiovascular_disease", "Autoimmune", "Acid_reflux", "Lung_disease") phenotype_summary <- data.frame( Condition = conditions, Recorded_cases = vapply( conditions, function(x) sum(meta[[x]] == 1, na.rm = TRUE), integer(1) ), Non_missing = vapply( conditions, function(x) sum(!is.na(meta[[x]])), integer(1) ) ) knitr::kable( phenotype_summary, caption = "Recorded phenotypes in the selected demonstration subset" )
This analysis has several important caveats:
The repository also contains the source for a separate 5,000-sample companion analysis. It is not installed as a package vignette.
The dCST framework provides a clustering-free approach to
community typing that is deterministic, hierarchical, and biologically
interpretable. The combination of depth-1 dCSTs (dominant species) and depth-2
refinement (co-dominance patterns) captures structure at multiple resolutions
without requiring any parameter tuning beyond the minimum support threshold
n0.
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.