Compositional data are first expressed in isometric logratio (pivot)
coordinates via pivotCoord. An alternating optimisation
scheme then jointly estimates cell-level contamination flags in the original
raw-part space, maps them to coordinate-level weights through the
contrast matrix, and computes a weighted PCA in the ilr space. Loadings
and the robust centre are back-transformed to the clr space for
interpretation, following the same convention as pcaCoDa.
By default the alternating scheme is initialised from a robust MCD fit of the
ilr coordinates (init = "mcd", the default in Templ 2026, underpinning
its breakdown results); a fast classical SVD start is available via
init = "classical". Robustness is controlled by a bounded loss
rho, which defaults to the Tukey biweight (95% Gaussian efficiency),
as in the paper.
The detection step uses a projection-based test derived from the
propagation theorem: if raw part \(j\) alone is contaminated, the
perturbation in ilr space lies along the direction
\(d_j = V^\top (e_j - 1/D)\) (normalised to unit length). MAD-standardised
ilr residuals are projected onto each \(d_j\), and the projection score is
standardised by its correlation-adjusted scale
\(\sigma_{p,j} = \sqrt{d_j^\top R\, d_j}\), where \(R\) is the correlation
matrix of the standardised residuals, so the standardised score is
approximately \(N(0,1)\) under the null. A cell is flagged when this score
exceeds the Bonferroni-adjusted threshold \(z_{1 - \alpha/(2D)}\).
The weight mapping from raw-part flags to ilr weights is fully vectorised
via a matrix multiply (no triple-nested loop). For each observation
\(i\) and ilr coordinate \(l\) the weight is
$$w_{il} = \max\Bigl(w_0,\;
\exp\bigl(\sum_j c_{ij}\,\log\bigl(1-(1-w_0)\,|V_{jl}|/\max_j |V_{jl}|\bigr)\bigr)\Bigr)$$
where \(c_{ij}\) are the binary contamination flags (1 = flagged).
The weighted covariance uses per-coordinate normalisation:
entry \((l, m)\) of the covariance is normalised by
\(\sqrt{\sum_i w_{il} \cdot \sum_i w_{im}}\) so that coordinates with
different total weights are treated correctly.
Implementation choices beyond Templ (2026). Two conveniences are
not part of the published algorithm and exist only for ease of use: (i) when
k is NULL it is chosen by the Kaiser criterion (retain
eigenvalues above the average) --- the paper leaves k to the analyst,
so set it explicitly (e.g. from a scree plot or a clear eigenvalue gap) for
publication-quality results; and (ii) the ridge parameter reg
(default 0, i.e. the paper's algorithm) adds a small diagonal term for
ill-conditioned high-dimensional compositions, a regime the paper lists as
future work. Convergence is declared when the objective's relative change
falls below tol, or when a period-2 limit cycle is detected (a few
borderline cells oscillating under the hard-threshold reweighting); in the
latter case it stops at the stabilised cycle (the solution is stable up to
those few cells) and cycled is TRUE.