cellNetCoDa: Cellwise-robust network estimation for compositional data

Description

Estimates a sparse conditional-independence (Gaussian graphical) network among the parts of a composition that is robust against cellwise contamination --- individual corrupted parts --- on top of the closure constraint of the simplex.

Usage

cellNetCoDa(
  X,
  method = c("mb", "glasso"),
  rho_grid = NULL,
  lambda_grid = NULL,
  B = 20,
  sub_ratio = 0.7,
  thresh = 0.05,
  k = NULL,
  rule = "and",
  taxa = colnames(X),
  ...
)

Value

An object of class "cellNetCoDa" with components:

adj: $D \times D$ symmetric 0/1 adjacency matrix with zero diagonal and dimnames given by taxa: the estimated conditional-independence graph of the parts
taxa: character vector of node (part) names
n_edges: number of edges (the upper-triangle sum of adj)
weights: $n \times D$ matrix of cellwise weights in $[0, 1]$ (1 = clean, smaller = more strongly down-weighted)
penalty: the StARS-selected penalty (rho for glasso, lambda for MB)
instability: the StARS instability curve over the penalty grid
grid: the penalty grid the StARS curve was evaluated on
thresh: the StARS instability threshold used for selection
method: the back end used, "mb" or "glasso"
flags: $n \times D$ logical matrix of cellwise contamination flags

Arguments

X: an $n \times D$ matrix or data.frame of compositional parts (columns are parts). Counts or proportions are accepted; count zeros are pre-treated by a Bayesian-multiplicative replacement (cmultRepl, method "GBM") and the rows closed to the simplex internally. At least three parts ($D \ge 3$).
method: back end for the sparse network estimate: "mb" (default) for cellwise-robust Meinshausen--Buhlmann neighbourhood selection in clr space, or "glasso" for the graphical lasso on the clr covariance.
rho_grid: numeric vector of graphical-lasso penalties (decreasing order recommended). Required when method = "glasso"; ignored otherwise.
lambda_grid: numeric vector of neighbourhood-selection penalties (decreasing). Required when method = "mb"; ignored otherwise.
B: number of subsamples used by the StARS penalty selection (default: 20).
sub_ratio: subsample fraction for StARS (default: 0.7).
thresh: StARS edge-instability threshold (default: 0.05). The densest graph whose monotonised instability stays at or below thresh is selected.
k: number of principal components passed to the cellwise detector cellPcaCoDa. If NULL (default) it is chosen there by the Kaiser criterion; set it explicitly (k < D - 1) for reproducible results.
rule: neighbourhood-combination rule for method = "mb": "and" (default; an edge requires both nodes to select each other) or "or". Ignored for "glasso".
taxa: optional character vector of node (part) names for the returned graph; defaults to colnames(X), or T1, ..., TD when X is unnamed.
...: further arguments for the cellwise detector, e.g. grading ("marginal", "projection" or "isotropic"), the flag floor w0, the detection cutoff k_thr, or the number of de-masking passes iter; remaining arguments are passed on to cellPcaCoDa. These are used both in the StARS subsampling and in the final fit.

Author

Matthias Templ

Details

The estimator chains four steps. (i) Count zeros are replaced by a Bayesian-multiplicative replacement and the rows are closed to the simplex. (ii) Cellwise outliers are flagged by a masking-resistant conditional-pivot test together with the projection flags of cellPcaCoDa and turned into per-cell weights; a cell-cleaned centred-logratio (clr) covariance is then formed by shrinking each flagged cell toward the robust clr centre, $$\tilde{C}_{ij} = w_{ij}\,C_{ij} + (1 - w_{ij})\,\mu_j,$$ where $C$ is the clr matrix, $w_{ij}\in[0,1]$ the cell weight and $\mu$ the robust clr centre. (iii) A sparse graph is read off the cleaned covariance, either by cellwise-robust Meinshausen--Buhlmann neighbourhood selection (method = "mb", the default) or by the graphical lasso on the clr covariance (method = "glasso"). (iv) The sparsity penalty is chosen by StARS.

Penalty selection (StARS). For each penalty on the supplied grid the graph is refitted on B random subsamples (fraction sub_ratio); the edge-selection frequency $\theta$ gives the edge instability $2\theta(1-\theta)$, averaged over edges. Following Liu et al. (2010) the instability is monotonised and the densest graph whose instability stays at or below thresh (the smallest stable penalty) is returned.

Back ends. "mb" regresses each part's clr value on the cell-cleaned clr of the remaining parts by a weighted lasso and combines the neighbourhoods by rule; "glasso" penalises the clr covariance directly and reads the graph from the exact zeros of the estimated precision matrix (the SPIEC-EASI convention). Both share the same subsampling, grid and instability convention, so penalties are comparable across back ends.

References

Liu, H., Roeder, K. and Wasserman, L. (2010). Stability approach to regularization selection (StARS) for high-dimensional graphical models. Advances in Neural Information Processing Systems 23.

Meinshausen, N. and Buhlmann, P. (2006). High-dimensional graphs and variable selection with the lasso. The Annals of Statistics 34(3), 1436--1462.

Friedman, J., Hastie, T. and Tibshirani, R. (2008). Sparse inverse covariance estimation with the graphical lasso. Biostatistics 9(3), 432--441.

Kurtz, Z. D., Muller, C. L., Miraldi, E. R., Littman, D. R., Blaser, M. J. and Bonneau, R. A. (2015). Sparse and compositionally robust inference of microbial ecological networks. PLoS Computational Biology 11(5), e1004226.

Raymaekers, J. and Rousseeuw, P. J. (2024). The cellwise minimum covariance determinant estimator. Journal of the American Statistical Association 119(548), 2610--2621.

Martin-Fernandez, J. A., Hron, K., Templ, M., Filzmoser, P. and Palarea-Albaladejo, J. (2015). Bayesian-multiplicative treatment of count zeros in compositional data sets. Statistical Modelling 15(2), 134--158.

Examples

Run this code

# \donttest{
data(expenditures)            # 5-part compositional budget data
set.seed(1)
## cellwise-robust network via MB neighbourhood selection (needs 'glmnet')
if (requireNamespace("glmnet", quietly = TRUE)) {
  net <- cellNetCoDa(expenditures, method = "mb",
                     lambda_grid = c(0.4, 0.3, 0.2, 0.1, 0.05),
                     B = 10, k = 2)
  net                         # print method: back end, edges, penalty
  net$adj                     # estimated 0/1 adjacency among the parts
}

## graphical-lasso back end, SPIEC-EASI style (needs 'glasso')
if (requireNamespace("glasso", quietly = TRUE)) {
  net_gl <- cellNetCoDa(expenditures, method = "glasso",
                        rho_grid = c(0.4, 0.3, 0.2, 0.1, 0.05),
                        B = 10, k = 2)
  net_gl$n_edges
}
# }

Run the code above in your browser using DataLab