API Reference

cca(X::Matrix{Float64}, Y::Matrix{Float64}; method::Symbol = :svd,
	outdim::Int = min(size(X,1), size(Y,1)))

Fit a canonical correlation analysis (CCA) between two sets of variables, after Weenink (2003)

Arguments

• X::Matrix{Float64}: 2d array of floats; the first variable set, with variables in ROWS and observations in COLUMNS (dx×n)
• Y::Matrix{Float64}: 2d array of floats; the second variable set (dy×n), with the same number of observations (columns) as X
• method::Symbol: The solver, :svd (works directly from the centered data, more stable) or :cov (works from covariance matrices). Both give the same result. Defaults to :svd
• outdim::Int: The number of canonical components to return; must be in 1:min(dx, dy). Defaults to min(dx, dy)

Value

A Cca holding the X and Y means, the canonical directions for each side, and the canonical correlations (descending). CCA finds pairs of directions, one in X-space and one in Y-space, whose projected variates are maximally correlated; the k-th pair is the most correlated subject to being uncorrelated with the previous k−1. The default solver follows the SVD-based approach of Weenink (2003), which avoids forming covariance matrices for numerical stability. Warns when n is not comfortably larger than dx and dy, where CCA becomes unstable

cca_transform(M::Cca, Z::AbstractMatrix, c::Symbol)

Project data onto the canonical directions of a fitted CCA model

Arguments

• M::Cca: A fitted CCA model, as returned by cca
• Z::AbstractMatrix: 2d array of floats; the data to project, with variables in rows and observations in columns, matching the side selected by c
• c::Symbol: Which side to project, :x (using the X directions and X mean) or :y (using the Y directions and Y mean)

Value

2d array of floats; the p×n matrix of canonical variates (each observation's coordinates along the canonical directions). Throws an ArgumentError if c is neither :x nor :y

jive(Xs::Vector{Matrix{Float64}}, r::Int, ri::Vector{Int}; kwargs...)

Positional-argument form of jive for supplying the joint rank r and individual ranks ri directly

Arguments

• Xs::Vector{Matrix{Float64}}: The data blocks, each pᵢ×n
• r::Int: The joint rank
• ri::Vector{Int}: The individual ranks, one per block
• kwargs...: Any keyword arguments accepted by the keyword form of jive

Value

A Jive; identical to calling jive(Xs; r = r, ri = ri, kwargs...)

jive(Xs::Vector{Matrix{Float64}}; r = nothing, ri = nothing, scale = true,
	 center = true, tol = nothing, maxiter = 1000, nperm = 100, alpha = 0.05)

Fit a JIVE (Joint and Individual Variation Explained) decomposition of several data blocks sharing a common set of observations

Arguments

• Xs::Vector{Matrix{Float64}}: The data blocks, each pᵢ×n with the SAME number of columns n (observations); rows (features) may differ per block
• r: The joint rank; if nothing, estimated by permutation test. Defaults to nothing
• ri: The vector of individual ranks; if nothing, estimated by permutation test. Defaults to nothing
• scale::Bool: Whether to scale each block by its Frobenius norm (and √(total elements)) so blocks contribute comparably. Defaults to true
• center::Bool: Whether to row-center each block. Defaults to true
• tol: The convergence threshold; if nothing, set to 1e-6·‖stacked data‖. Defaults to nothing
• maxiter::Int: The maximum core-fit iterations. Defaults to 1000
• nperm::Int: The permutations per rank test (when ranks are estimated). Defaults to 100
• alpha::Float64: The significance level for the rank test. Defaults to 0.05

Value

A Jive with the joint and individual decompositions. Each block is first row-centered and Frobenius-scaled so no block dominates by sheer magnitude. If the joint and individual ranks are not supplied, they are estimated by a permutation test; the decomposition is then computed by alternating between estimating the shared joint structure and each block's individual structure until they converge. The result splits each block into joint structure (a common row space across all blocks) and block-specific individual structure orthogonal to it, after Lock et al. (2013)

nmf(X::Matrix{Float64}; k::Int = minimum(size(X)), init::Symbol = :auto,
	tol::Real = 1e-4, maxiter::Int = 200, alpha_w::Real = 0.0,
	alpha_h::Union{Nothing,Real} = nothing, l1_ratio::Real = 0.0,
	shuffle::Bool = false, rng::AbstractRNG = Random.GLOBAL_RNG,
	w_init = nothing, h_init = nothing)

Fit a nonnegative matrix factorization using Fast HALS coordinate descent

Arguments

• X::Matrix{Float64}: 2d array of nonnegative finite floats; observations in rows and features in columns
• k::Int: The factorization rank. Defaults to min(n,p)
• init::Symbol: Initialization method: :auto, :random, :nndsvd, :nndsvda, :nndsvdar, or :custom. :auto uses :nndsvda when possible and :random otherwise. Defaults to :auto
• tol::Real: Stopping tolerance for the normalized projected-gradient violation. Defaults to 1e-4
• maxiter::Int: Maximum number of alternating coordinate-descent iterations. Defaults to 200
• alpha_w::Real: Elastic-net strength for W, scaled internally by the number of features. Defaults to 0
• alpha_h::Union{Nothing,Real}: Elastic-net strength for H, scaled internally by the number of observations. nothing uses alpha_w. Defaults to nothing
• l1_ratio::Real: Fraction of regularization assigned to L1, in [0,1]; the remainder is squared Frobenius regularization. Defaults to 0
• shuffle::Bool: Whether to randomize the component-coordinate order at every W and H update. Defaults to false
• rng::AbstractRNG: Random-number generator used by randomized initialization and coordinate shuffling. Defaults to Julia's global RNG
• w_init: Custom n×k nonnegative initial W, required only with init=:custom
• h_init: Custom k×p nonnegative initial H, required only with init=:custom

Value

An Nmf holding nonnegative factors W and H, their Frobenius reconstruction error, convergence information, and regularization settings. It minimizes the scikit-learn NMF objective: one-half squared reconstruction error plus size-scaled elastic-net penalties on W and H, using alternating Fast HALS coordinate descent. Warns if maxiter is reached before convergence

nmf_invtransform(M::Nmf, W::Matrix{Float64})

Reconstruct observations in the original feature space from nonnegative NMF scores

Arguments

• M::Nmf: A fitted NMF model, as returned by nmf
• W::Matrix{Float64}: 2d array of floats; the n×k score matrix to reconstruct, with k matching the number of fitted components

Value

2d array of floats; the n×p reconstruction W*M.h

nmf_transform(M::Nmf, X::Matrix{Float64}; tol::Real = 1e-4,
	maxiter::Int = 200, shuffle::Bool = false,
	rng::AbstractRNG = Random.GLOBAL_RNG)

Represent new nonnegative observations using the fixed components of a fitted NMF

Arguments

• M::Nmf: A fitted NMF model, as returned by nmf
• X::Matrix{Float64}: 2d array of nonnegative finite floats; new observations in rows and the same p features used to fit M in columns
• tol::Real: The stopping tolerance for the normalized projected-gradient violation. Defaults to 1e-4
• maxiter::Int: The maximum number of coordinate-descent sweeps. Defaults to 200
• shuffle::Bool: Whether to randomize the component-update order in each sweep. Defaults to false
• rng::AbstractRNG: Random-number generator used only when shuffle=true. Defaults to Julia's global RNG

Value

2d array of floats; the n×k nonnegative score matrix W minimizing ‖X − W*M.h‖² with the fitted components held fixed and the model's W regularization applied

pca(X::Matrix{Float64}; k::Int = minimum(size(X)),
	standardize::Bool = false, method::Symbol = :auto)

Fit a principal component analysis (PCA) model

Arguments

• X::Matrix{Float64}: 2d array of floats; observations (rows) by features (columns)
• k::Int: The number of components to retain; 1 ≤ k ≤ min(n,p). Defaults to min(n,p)
• standardize::Bool: Whether to scale columns to unit variance in addition to centering. Defaults to false (center only)
• method::Symbol: :auto (pick by shape — :cov when n ≥ p, :svd when p > n), :cov (eigendecomposition of the p×p covariance), or :svd (SVD of the centered data). Defaults to :auto

Value

A Pca holding the column means, scales, k loadings, component variances, and proportion of variance explained

pca_invtransform(m::Pca, scores::Matrix{Float64})

Reconstruct data in the original feature space from principal-component scores

Arguments

• m::Pca: A fitted PCA model, as returned by pca
• scores::Matrix{Float64}: 2d array of floats; the n×k matrix of principal-component scores to invert, with k matching the number of components retained in m

Value

2d array of floats; the n×p reconstruction in the original units. Exact only when all components are retained (k = p); otherwise a low-rank approximation

pca_transform(m::Pca, X::Matrix{Float64})

Project data onto the principal directions of a fitted PCA model

Arguments

• m::Pca: A fitted PCA model, as returned by pca
• X::Matrix{Float64}: 2d array of floats; the observations (rows) by features (columns) to project, with the same features as the training data

Value

2d array of floats; the n×k matrix of principal-component scores

plsda(X::Matrix{Float64}, y::Vector, ncomp::Int;
	  scale = true, tol = 1e-6, max_iter = 100, levels = nothing)

Fit a PLS discriminant analysis (PLS-DA) for multiclass problems — the non-sparse counterpart of splsda. The class labels are dummy-encoded and the problem solved as a PLS regression of X onto the indicator matrix: each component is the rank-1 approximation of the cross-product Mₕ = XₕᵀYₕ by alternating power iteration, after which both X and Y are deflated by regression on the X-variate (mixOmics' regression mode). Unlike splsda, no L1 penalty is applied, so every variable contributes to each loading.

plskern(X::Matrix{Float64}, Y::Matrix{Float64}; nlv::Int = 2,
		standardize::Bool = false, method::Symbol = :algo1)

Fit a kernel partial least squares (PLS) regression using the improved algorithms of Dayal & MacGregor (1997). The caller's X and Y are left untouched; centering and scaling are done on internal copies

Arguments

• X::Matrix{Float64}: 2d array of floats; the n×p predictor matrix
• Y::Matrix{Float64}: 2d array of floats; the n×q response matrix
• nlv::Int: The number of latent variables; clamped to min(nlv, n, p). Defaults to 2
• standardize::Bool: Whether to scale each column to unit standard deviation in addition to centering. Defaults to false (center only)
• method::Symbol: :algo1 (kernel via XᵀY) or :algo2 (kernel via XᵀX, faster when p ≪ n). Both give identical results. Defaults to :algo1

Value

A Plskern holding the weight (W, R), loading (P, Q), and score (T) matrices together with the centering means and scales

plskern_coef(m::Plskern; nlv::Int = size(m.R, 2))

Assemble the regression coefficient matrix and intercept from a fitted PLS model

Arguments

• m::Plskern: A fitted PLS model, as returned by plskern
• nlv::Int: The number of latent variables to include; clamped to the number fitted. Defaults to all fitted components

Value

A tuple (B, intercept) where B is the p×q coefficient matrix and intercept is the 1×q intercept, such that a prediction is intercept .+ Xnew * B. The scales are folded in so B maps raw (uncentered, unscaled) X to Y

plskern_predict(m::Plskern, Xnew::Matrix{Float64}; nlv::Int = size(m.R, 2))

Predict responses for new observations from a fitted PLS model

Arguments

• m::Plskern: A fitted PLS model, as returned by plskern
• Xnew::Matrix{Float64}: 2d array of floats; the new observations (rows) by predictors (columns), with the same p predictors as the training data
• nlv::Int: The number of latent variables to use; clamped to the number fitted. Defaults to all fitted components

Value

2d array of floats; the predicted n×q responses, intercept .+ Xnew * B

plskern_transform(m::Plskern, Xnew::Matrix{Float64}; nlv::Int = size(m.R, 2))

Project new observations onto the PLS latent space (compute their X-scores)

Arguments

• m::Plskern: A fitted PLS model, as returned by plskern
• Xnew::Matrix{Float64}: 2d array of floats; the new observations (rows) by predictors (columns), with the same p predictors as the training data
• nlv::Int: The number of latent variables to project onto; clamped to the number fitted. Defaults to all fitted components

Value

2d array of floats; the n×nlv matrix of X-scores. For the training data this reproduces the stored scores m.T, since scores are linear: T = Xc·R

pmd(X::Matrix{Float64}; sumabs::Real = 0.4,
	sumabsu::Union{Nothing,Real} = nothing,
	sumabsv::Union{Nothing,Real} = nothing,
	K::Int = 1, niter::Int = 20, center::Bool = true)

Fit a penalized matrix decomposition with L1 penalties on both factors, PMD(L1,L1)

Arguments

• X::Matrix{Float64}: 2d array of floats; the n×p data matrix
• sumabs::Real: The shared sparsity budget in (0,1], scaled internally to give the per-column L1 budgets su = √n·sumabs on u and sv = √p·sumabs on v. Smaller values give sparser factors. Ignored when both sumabsu and sumabsv are supplied. Defaults to 0.4
• sumabsu::Union{Nothing,Real}: An explicit L1 budget on each u column, in [1, √n]; overrides sumabs when given together with sumabsv. Defaults to nothing
• sumabsv::Union{Nothing,Real}: An explicit L1 budget on each v column, in [1, √p]; overrides sumabs when given together with sumabsu. Defaults to nothing
• K::Int: The number of components to extract; must be in 1:min(n,p). Defaults to 1
• niter::Int: The maximum number of power-iteration steps per component. Defaults to 20
• center::Bool: Whether to subtract the grand mean of X before decomposing. Defaults to true

Value

A Pmd holding the n×K left factors u, the p×K right factors v, the K weights d, the resolved budgets su and sv, the number of components K, and the grand mean removed (or NaN when center=false)

scca(X::Matrix{Float64}, Y::Matrix{Float64}; penaltyx::Real = 0.3,
	 penaltyz::Real = 0.3, K::Int = 1, niter::Int = 15, standardize::Bool = true)

Fit a sparse canonical correlation analysis (sparse CCA) between two sets of variables, after Witten, Tibshirani & Hastie (2009)

Arguments

• X::Matrix{Float64}: 2d array of floats; the first variable set, with variables in ROWS and observations in COLUMNS (dx×n)
• Y::Matrix{Float64}: 2d array of floats; the second variable set (dy×n), same number of observations (columns) as X
• penaltyx::Real: The L1 penalty on the X canonical vectors, in (0,1]; smaller values give sparser vectors. Defaults to 0.3
• penaltyz::Real: The L1 penalty on the Y canonical vectors, in (0,1]. Defaults to 0.3
• K::Int: The number of canonical components to extract; must be in 1:min(dx,dy). Defaults to 1
• niter::Int: The maximum power-iteration steps per component. Defaults to 15
• standardize::Bool: Whether to center and scale each variable to unit variance before fitting. Defaults to true

Value

An Scca holding the sparse X and Y canonical vectors (u, v), the component weights d, the sample correlations of the paired variates, the penalties, and K. Sparse CCA finds pairs of sparse vectors, one in X-space and one in Y-space, whose projected variates are maximally correlated under an L1 penalty (so each vector selects only a few variables). It uses the diagonal-covariance penalized formulation of Witten, Tibshirani & Hastie (2009): each component is a rank-1 penalized approximation of XᵀY solved by alternating soft-thresholded power iteration, after which the data are deflated before the next component

spc(X::Matrix{Float64}; k::Int = 2, c::Real = sqrt(size(X,2))/2,
	standardize::Bool = false, tol::Real = 1e-7, niter::Int = 20)

Fit a sparse principal component analysis using the SPC criterion (deflation variant)

Arguments

• X::Matrix{Float64}: 2d array of floats; the n×p data matrix
• k::Int: The number of sparse components to extract. Defaults to 2
• c::Real: The L1 budget on each loading, in [1, √p]; smaller values give sparser loadings, c = √p gives ordinary (dense) PCA. Defaults to √p / 2
• standardize::Bool: Whether to scale each column to unit standard deviation in addition to centering. Defaults to false (center only)
• tol::Real: The ℓ₁ convergence tolerance per component. Defaults to 1e-7
• niter::Int: The maximum number of power-iteration steps per component. Defaults to 20

Value

An Spc holding the column means, scales, k sparse loadings, component variances, and cumulative proportion of variance explained. Components are found sequentially: each loading is fit on the residual, then removed by rank-1 deflation before the next is extracted

spc_orth(X::Matrix{Float64}; k::Int = 2, c::Real = sqrt(size(X,2))/2,
		 standardize::Bool = false, tol::Real = 1e-6, niter::Int = 20)

Fit a sparse principal component analysis with orthogonal scores (orthogonal variant of the SPC criterion)

Arguments

• X::Matrix{Float64}: 2d array of floats; the n×p data matrix
• k::Int: The number of sparse components to extract. Defaults to 2
• c::Real: The L1 budget on each loading, in [1, √p]; smaller values give sparser loadings. Defaults to √p / 2
• standardize::Bool: Whether to scale each column to unit standard deviation in addition to centering. Defaults to false (center only)
• tol::Real: The ℓ₁ convergence tolerance per component. Defaults to 1e-6
• niter::Int: The maximum number of power-iteration steps per component. Defaults to 20

Value

An Spc holding the column means, scales, k sparse loadings, component variances, and cumulative proportion of variance explained. Unlike spc, each component's scores are constrained to be orthogonal to all previous components' scores (via projection rather than deflation), so the score directions form an orthonormal set

splsda(X::Matrix{Float64}, y::Vector, ncomp::Int, keepX::Vector{Int};
	   scale = true, tol = 1e-6, max_iter = 100, levels = nothing)

Fit a sparse PLS discriminant analysis (sPLS-DA) for multiclass problems, after Lê Cao, Boitard & Besse (2011)

Arguments

• X::Matrix{Float64}: 2d array of floats; the n×p predictor matrix
• y::Vector: The class labels, one per observation (n total)
• ncomp::Int: The number of components to extract
• keepX::Vector{Int}: The number of variables to keep per component; must have length ncomp (one sparsity budget per component)
• scale::Bool: Whether to scale columns to unit variance in addition to centering. Defaults to true
• tol::Real: The convergence tolerance on the loading iterates. Defaults to 1e-6
• max_iter::Int: The maximum number of inner iterations per component. Defaults to 100
• levels: An optional class ordering passed to the dummy encoding; fixes the column order of the classes. Defaults to nothing (sorted unique labels)

Value

An Splsda holding the X and Y variates (scores), the sparse X and dense Y loadings, the number of components, the per-component keepX, the one-hot class encoding, and the class labels. The class labels are dummy-encoded into an indicator matrix and the problem solved as a sparse PLS regression of X onto that matrix: each component is the penalized rank-1 approximation of the cross-product Mₕ = XₕᵀYₕ, solved by alternating soft-thresholded power iteration, after which X (only) is deflated by regression on its variate (Y is not deflated — the discriminant regression mode)

vip(m)

Compute Variable Importance in Projection (VIP) scores for a fitted PLS discriminant model

Arguments

• m: a fitted model (Plsda or Splsda) exposing the X-loadings (loadings_X), the X-variates (variates_X), the dummy-encoded response (Y_dummy), and the number of components (ncomp)

Value

A p×ncomp matrix of VIP scores, where p is the number of variables (rows of the X-loadings). Entry [j, h] is the importance of variable j using the first h components; the last column is the cross-component VIP that summarizes each variable's contribution over all components. The scores are normalized so that the mean squared VIP equals 1, so a value above 1 marks an above-average-importance variable — the conventional selection threshold. This follows the VIP definition of Wold et al. as implemented in mixOmics::vip: each component's contribution is weighted by the redundancy (the squared correlation between the dummy response and that component's variate) it explains, so components that better separate the classes count more toward the score

Cca{T}

Container for a fitted canonical correlation analysis, as returned by cca

Fields

• xmean::Vector{T}: The mean of each X variable (length dx), removed during centering
• ymean::Vector{T}: The mean of each Y variable (length dy), removed during centering
• xproj::Matrix{T}: The dx×p X canonical directions (columns), one per component
• yproj::Matrix{T}: The dy×p Y canonical directions (columns), one per component
• corrs::Vector{T}: The p canonical correlations, in [0,1] and descending
• nobs::Int: The number of observations used (or −1 when fit via the covariance path, which doesn't carry n through)

Jive{T}

Container for a fitted JIVE (Joint and Individual Variation Explained) decomposition, as returned by jive

Fields

• J::Vector{Matrix{T}}: The joint structure, one pᵢ×n matrix per data block; all blocks share the same n-dimensional joint row space
• A::Vector{Matrix{T}}: The individual structure, one pᵢ×n matrix per block, orthogonal to the joint structure
• S::Matrix{T}: The r×n joint scores (orthonormal rows), the basis of the shared joint row space
• U::Vector{Matrix{T}}: The joint loadings, one pᵢ×r block, so that J[i] = U[i]·S
• Si::Vector{Matrix{T}}: The individual scores, one rᵢ×n matrix per block
• Wi::Vector{Matrix{T}}: The individual loadings, one pᵢ×rᵢ matrix per block, so that A[i] = Wi[i]·Si[i]
• r::Int: The joint rank (dimension of the shared structure)
• ri::Vector{Int}: The individual ranks, one per block

Nmf{T}

Container for a fitted nonnegative matrix factorization, as returned by nmf

Fields

• w::Matrix{T}: The n×k nonnegative sample scores, one row per observation
• h::Matrix{T}: The k×p nonnegative components, one column per input feature
• reconstruction_error::T: The Frobenius norm ‖X − w*h‖_F of the final reconstruction residual
• niter::Int: The number of alternating coordinate-descent iterations performed
• converged::Bool: Whether the normalized projected-gradient violation reached the requested tolerance before maxiter
• alpha_w::T: The elastic-net regularization strength used for w
• alpha_h::T: The elastic-net regularization strength used for h
• l1_ratio::T: The fraction of regularization assigned to the L1 penalty, in [0,1]; the remainder is assigned to the squared Frobenius penalty

Pca{T}

Container for a fitted PCA model, as returned by pca

Fields

• mean::Vector{T}: The column means of the training data (length p), removed during centering
• scale::Vector{T}: The column scales (length p) — the column standard deviations when standardize=true, otherwise ones
• loadings::Matrix{T}: The p×k principal directions (right singular vectors of the centered data / leading eigenvectors of the covariance matrix), unit-norm columns ordered by decreasing explained variance
• variances::Vector{T}: The variance explained by each of the k components (the eigenvalues, σᵢ² / (n-1))
• propOFvar::Vector{T}: The proportion of total variance explained by each component (variances ./ total variance)

Plsda{T}

Container for a fitted PLS discriminant analysis, as returned by plsda.

Fields

• variates_X::Matrix{T}: n×ncomp X-scores (sample coordinates)
• variates_Y::Matrix{T}: n×ncomp Y-scores
• loadings_X::Matrix{T}: p×ncomp dense X-loadings
• loadings_Y::Matrix{T}: K×ncomp Y-loadings (K = number of classes)
• ncomp::Int: number of components
• Y_dummy::Matrix{T}: n×K one-hot encoding of the labels
• classes::Vector: class labels, in Y_dummy column order

Plskern{T}

Container for a fitted kernel PLS regression model, as returned by plskern. Holds the standard PLS factor matrices and the centering/scaling statistics

Fields

• W::Matrix{T}: The p×nlv weight matrix (unit-norm columns)
• P::Matrix{T}: The p×nlv X-loadings
• Q::Matrix{T}: The q×nlv Y-loadings
• R::Matrix{T}: The p×nlv weight matrix in the deflated basis, so that scores satisfy T = Xc·R (the SIMPLS-style weights)
• T::Matrix{T}: The n×nlv X-scores (orthogonal columns)
• xmeans::Vector{T}: The column means of X removed during centering
• xscales::Vector{T}: The column scales of X (standard deviations if standardize=true, otherwise ones)
• ymeans::Vector{T}: The column means of Y removed during centering
• yscales::Vector{T}: The column scales of Y

Pmd{T}

Container for a fitted penalized matrix decomposition, as returned by pmd

Fields

• u::Matrix{T}: The n×K left factors (one unit-ℓ₂, L1-penalized column per component)
• v::Matrix{T}: The p×K right factors (one unit-ℓ₂, L1-penalized column per component)
• d::Vector{T}: The K nonnegative singular-value-like weights, dₖ = uₖᵀ X vₖ
• sumabsu::Float64: The resolved L1 budget on each u column (in [1, √n])
• sumabsv::Float64: The resolved L1 budget on each v column (in [1, √p])
• K::Int: The number of components
• meanx::Float64: The grand mean removed during centering, or NaN if center=false

Scca{T}

Container for a fitted sparse canonical correlation analysis, as returned by scca

Fields

• u::Matrix{T}: The dx×K sparse X canonical vectors (one unit-ℓ₂, L1-penalized column per component); nonzero entries are the selected X features
• v::Matrix{T}: The dy×K sparse Y canonical vectors; nonzero entries are the selected Y features
• d::Vector{T}: The K weights dₖ = uₖᵀ Xᵀ Z vₖ (the penalized-covariance objective values)
• cors::Vector{T}: The K sample correlations between the paired canonical variates
• penaltyx::Float64: The L1 penalty applied to the X vectors, in (0,1]
• penaltyz::Float64: The L1 penalty applied to the Y vectors, in (0,1]
• K::Int: The number of components

Spc{T}

Container for a fitted sparse principal component analysis, as returned by spc or spc_orth

Fields

• mean::Vector{T}: The column means of the training data (length p), removed during centering
• scale::Vector{T}: The column scales (length p) — the column standard deviations when standardize=true, otherwise ones
• loadings::Matrix{T}: The p×k sparse loadings (one unit-ℓ₂, L1-penalized column per component), signs fixed by _sign_consistency_opt!
• variances::Vector{T}: The variance associated with each of the k components (dₖ² / (n-1), where dₖ is the component weight)
• propOFvar::Vector{T}: The cumulative proportion of variance explained by the first k components (see _prop_var_explained)

Splsda{T}

Container for a fitted sparse PLS discriminant analysis, as returned by splsda

Fields

• variates_X::Matrix{T}: The n×ncomp X-scores (sample coordinates for plotting/classification)
• variates_Y::Matrix{T}: The n×ncomp Y-scores (the class-side variates)
• loadings_X::Matrix{T}: The p×ncomp sparse X-loadings; the nonzero entries per column are the variables selected for that component
• loadings_Y::Matrix{T}: The K×ncomp Y-loadings (K = number of classes)
• ncomp::Int: The number of components fitted
• keepX::Vector{Int}: The number of variables kept per component (the sparsity budget)
• Y_dummy::Matrix{T}: The n×K one-hot encoding of the class labels
• classes::Vector: The class labels, in the column order of Y_dummy

_binary_search_opt(argu, sumabs)

Find the soft-threshold that makes a soft-thresholded vector meet a target L1 budget

Arguments

• argu: 1d array of floats; the vector to be thresholded
• sumabs: Float; the target ℓ₁/ℓ₂ ratio (the L1 budget)

Value

Float; the threshold d ≥ 0 such that ‖S(argu,d)‖₁/‖S(argu,d)‖₂ ≈ sumabs. Returns 0 when the vector already satisfies the budget. Implements the bisection step bracketing d in [0, max|argu|] and halving for up to 150 iterations or until the bracket is below 1e-6

_cca_cov_opt(Cxx, Cyy, Cxy, xmean, ymean, p::Int)

Solve CCA from the covariance matrices via a generalized eigenproblem

Arguments

• Cxx: 2d array of floats; the dx×dx X covariance
• Cyy: 2d array of floats; the dy×dy Y covariance
• Cxy: 2d array of floats; the dx×dy cross-covariance
• xmean: 1d array of floats; the X means, stored in the result
• ymean: 1d array of floats; the Y means, stored in the result
• p::Int: The number of canonical components to return

Value

A Cca with the top-p canonical directions and correlations (its nobs is −1, since this path works from covariances and doesn't carry n). Solves the classical canonical-correlation generalized eigenproblem (Weenink 2003), reducing through whichever side is smaller (dx ≤ dy or dx > dy) for efficiency: the squared canonical correlations are the eigenvalues, one side's directions are the eigenvectors, and the other side's are recovered and normalized to unit variance by _qnormalize!

_cca_svd_opt(Zx, Zy, xmean, ymean, p::Int)

Solve CCA directly from the centered data via SVD (no covariance matrices formed)

Arguments

• Zx: 2d array of floats; the centered X data (dx×n), OVERWRITTEN by svd!
• Zy: 2d array of floats; the centered Y data (dy×n), OVERWRITTEN by svd!
• xmean: 1d array of floats; the X means, stored in the result
• ymean: 1d array of floats; the Y means, stored in the result
• p::Int: The number of canonical components to return

Value

A Cca with the top-p canonical directions and correlations. Works from the SVDs of Zx and Zy directly rather than forming covariance matrices, which is the numerically stable route advocated by Weenink (2003): the canonical correlations are the singular values of the overlap between the two right-singular-vector bases, and the canonical directions are recovered by rescaling the left singular vectors so the variates have unit variance

_center_scale(M::AbstractMatrix; scale = true)

Center, and optionally scale, the columns of a matrix

Arguments

• M::AbstractMatrix: 2d array of floats; the matrix to center/scale
• scale::Bool: Whether to divide each column by its (n−1) standard deviation after centering. Defaults to true

Value

2d array of floats; the column-centered (and, if scale=true, unit-variance) matrix. Constant columns (zero standard deviation) are guarded: their scale is treated as 1 and the column is set to all zeros rather than producing NaN. The n−1 (corrected) standard deviation matches mixOmics' colSds

_fast_init_v(x, z, K)

Initialize the Y canonical vectors for wide data via an SVD construction

Arguments

• x: 2d array of floats; the (observations × features) X block
• z: 2d array of floats; the (observations × features) Z block
• K: Int; the number of components to initialize

Value

2d array of floats; the K leading left singular vectors of zᵀ·(xxᵀ)^½, used as the starting v for the sparse iteration when both blocks are wider than they are tall. Building the init from xxᵀ (observations × observations) keeps the SVD small in the high-dimensional case

_finding_v!(v, s, z, c)

Project a vector onto the intersection of the unit ℓ₂ ball and the L1 budget, writing the result in place

Arguments

• v: 1d array of floats; the destination buffer for the result, overwritten in place
• s: 1d array of floats; preallocated scratch for the soft-thresholded vector
• z: 1d array of floats; the raw direction to project (e.g. Xᵀu)
• c: Float; the L1 budget (in [1, √p]); the returned v satisfies ‖v‖₂ = 1 and ‖v‖₁ ≤ c

Value

The buffer v, holding the unit-ℓ₂ vector closest to z whose ℓ₁ norm meets the budget. If the plain unit-normalized z already satisfies the budget it is returned unchanged; otherwise a soft-threshold δ is found by bisection (up to 150 steps) so that the normalized, thresholded vector has ‖v‖₁ ≈ c. This is the penalized update of the SPC criterion in Witten, Tibshirani & Hastie (2009)

_init_rsv(Xc, k::Int)

Compute the top-k right singular vectors of a centered data matrix, used to initialize the sparse power iteration

Arguments

• Xc: 2d array of floats; the centered (and optionally scaled) data matrix
• k::Int: The number of leading right singular vectors to return

Value

2d array of floats; the p×k matrix of leading right singular vectors of Xc. Two branches keep the eigenproblem small: for tall data (p ≤ n) it eigendecomposes the p×p Gram matrix XcᵀXc directly, otherwise it eigendecomposes the n×n Gram matrix XcXcᵀ and back-projects. Both use a partial eigendecomposition (only the top k eigenpairs) and return columns in descending eigenvalue order

_jive_perm_ranks_opt(Xc::Vector{Matrix{Float64}}, n::Int; nperm::Int,
					 alpha::Float64, conv::Float64, maxiter::Int,
					 maxrounds::Int = 10)

Estimate the joint and individual ranks by a permutation test

Arguments

• Xc::Vector{Matrix{Float64}}: The preprocessed data blocks, each pᵢ×n
• n::Int: The shared number of columns (observations)
• nperm::Int: The number of permutations used to build each null distribution
• alpha::Float64: The significance level; a singular value counts as signal if it exceeds the (1−alpha) quantile of the permuted null
• conv::Float64: The convergence threshold passed to the core fit
• maxiter::Int: The maximum core-fit iterations
• maxrounds::Int: The maximum number of estimate-then-refit rounds. Defaults to 10

Value

A tuple (rJ, rA) of the estimated joint rank and the vector of individual ranks. The joint rank is the number of leading singular values of the stacked (data − individual) residual that exceed the column-permuted null; each individual rank is found the same way on each block's (data − joint) residual. Ranks and the fit are alternated until the rank estimates stabilize, following the r.jive permutation procedure of Lock et al. (2013)

_jive_rjive_core_opt2(Xc::Vector{Matrix{Float64}}, n::Int, r::Int,
					  ri::Vector{Int}; conv::Float64, maxiter::Int)

Compute the JIVE decomposition for GIVEN joint and individual ranks, by alternating between estimating the joint and individual structure

Arguments

• Xc::Vector{Matrix{Float64}}: The preprocessed (centered/scaled) data blocks, each pᵢ×n with a shared column dimension n
• n::Int: The number of columns (observations) shared across blocks
• r::Int: The joint rank
• ri::Vector{Int}: The individual ranks, one per block
• conv::Float64: The convergence threshold on the change in J and A between iterations
• maxiter::Int: The maximum number of outer iterations

Value

A Jive with the joint (J, S, U) and individual (A, Si, Wi) parts. Implements the r.jive "orthIndiv" algorithm of Lock et al. (2013): wide blocks are first compressed to their score space (via SVD) for efficiency, then the joint structure is the rank-r SVD of the stacked residual (data minus individual), and each block's individual structure is the rank-rᵢ SVD of its residual projected off the joint row space and off the other blocks' individual spaces. The two estimates are alternated to convergence, then mapped back to the original row space

_l1_diff(a, b)

Compute the ℓ₁ distance between two vectors

Arguments

• a: 1d array of floats
• b: 1d array of floats

Value

Float; ‖a − b‖₁. Used as the convergence criterion between successive iterates

_l1_of_norm(a)

Compute the ℓ₁/ℓ₂ ratio of a vector

Arguments

• a: 1d array of floats

Value

Float; ‖a‖₁ / ‖a‖₂. This ratio is the quantity the L1 budget constrains; it ranges from 1 (maximally sparse) to √length(a) (fully dense)

_l1_of_norm_soft(a, d)

Compute the ℓ₁/ℓ₂ ratio of a vector AFTER soft-thresholding it by d, without materializing the thresholded vector

Arguments

• a: 1d array of floats; the unthresholded vector
• d: Float; the soft-threshold to apply

Value

Float; ‖S(a,d)‖₁ / ‖S(a,d)‖₂ where S is the soft-threshold operator. Used inside the binary search to test whether a candidate threshold meets the L1 budget

_l1diff(a, b)

Compute the ℓ₁ distance between two vectors

Arguments

• a: 1d array of floats
• b: 1d array of floats

Value

Float; ‖a − b‖₁. Used as the convergence criterion between successive iterates

_l2n_val(a)

Compute the ℓ₂ norm of a vector with a zero-guard

Arguments

• a: 1d array of floats

Value

Float; the Euclidean norm ‖a‖₂, or 0.05 when a is all zeros. The guard mirrors PMA's convention and prevents division by zero in the normalization steps

_matsqrt(A)

Compute the symmetric square root of a symmetric positive-semidefinite matrix

Arguments

• A: 2d array of floats; a symmetric (positive-semidefinite) matrix

Value

2d array of floats; the symmetric matrix R such that R·R = A, formed from the eigendecomposition with eigenvalues floored at zero (so tiny negative roundoff eigenvalues don't produce complex results)

_nmf_cd_sweep!(A, F, G, gram, cross, gradient, order,
	l1_regularization::Float64, l2_regularization::Float64)

Perform one Fast HALS coordinate-descent sweep over a nonnegative factor

Arguments

• A: 2d array of floats; the m×q target matrix
• F: 2d array of floats; the mutable m×k factor updated in place
• G: 2d array of floats; the fixed q×k factor for this half-step
• gram: Preallocated k×k buffer for GᵀG
• cross: Preallocated m×k buffer for A*G
• gradient: Preallocated length-m gradient buffer
• order: Length-k integer vector giving the coordinate-update order
• l1_regularization::Float64: L1 term added to every coordinate gradient
• l2_regularization::Float64: L2 term added to the Gram-matrix diagonal

Value

Float; the sum of absolute projected-gradient violations before each coordinate update. The sweep updates every column of F in place using max(F[:,t] - gradient/GᵀG[t,t], 0). All matrix products and work vectors are reused, so the iterative hot path does not allocate an m×q reconstruction

_nmf_fit_cd!(X, w, ht, l1_w::Float64, l1_h::Float64,
	l2_w::Float64, l2_h::Float64; tol::Float64, maxiter::Int,
	shuffle::Bool, rng::AbstractRNG, update_h::Bool = true)

Fit mutable NMF factors with alternating Fast HALS coordinate descent

Arguments

• X: 2d array of floats; the n×p nonnegative target matrix
• w: Mutable n×k score matrix
• ht: Mutable p×k transposed-component matrix; held fixed when update_h=false
• l1_w::Float64: Scaled L1 regularization applied to W
• l1_h::Float64: Scaled L1 regularization applied to H
• l2_w::Float64: Scaled L2 regularization applied to W
• l2_h::Float64: Scaled L2 regularization applied to H
• tol::Float64: Stopping tolerance for violation/initial_violation
• maxiter::Int: Maximum number of alternating sweeps
• shuffle::Bool: Whether to shuffle component coordinates before each half-step
• rng::AbstractRNG: Random-number generator used when shuffling
• update_h::Bool: Whether to update H as well as W. Defaults to true

Value

A tuple (niter, converged). The convergence rule matches scikit-learn's coordinate-descent NMF: stop when the projected-gradient violation relative to the first iteration is no greater than tol

_nmf_initialize(X::Matrix{Float64}, k::Int, init::Symbol,
	rng::AbstractRNG, w_init, h_init)

Resolve and construct the requested NMF initialization

Arguments

• X::Matrix{Float64}: 2d array of nonnegative floats; the n×p matrix to factorize
• k::Int: The requested factorization rank
• init::Symbol: :auto, :random, :nndsvd, :nndsvda, :nndsvdar, or :custom
• rng::AbstractRNG: Random-number generator for randomized initializations
• w_init: Optional n×k custom initial score matrix
• h_init: Optional k×p custom initial component matrix

Value

A tuple (w, ht) containing mutable Float64 factors. :auto selects :nndsvda when k <= min(n,p), otherwise :random. Throws an ArgumentError or DimensionMismatch for an invalid initialization request

_nmf_nndsvd_init(X::Matrix{Float64}, k::Int, variant::Symbol,
	rng::AbstractRNG; eps::Float64 = 1e-6)

Initialize NMF factors with Nonnegative Double Singular Value Decomposition

Arguments

• X::Matrix{Float64}: 2d array of nonnegative floats; the n×p matrix to factorize
• k::Int: The requested rank, no larger than min(n,p)
• variant::Symbol: :nndsvd to retain zeros, :nndsvda to replace zeros by mean(X), or :nndsvdar to replace zeros by small random values
• rng::AbstractRNG: Random-number generator used only by :nndsvdar
• eps::Float64: Values below this threshold are set to zero before applying the selected zero-filling rule. Defaults to 1e-6

Value

A tuple (w, ht) containing the n×k scores and p×k transposed components. Positive and negative singular-vector parts are compared without materializing four temporary vectors, and only the selected pair is written into the factors

_nmf_random_init(X::Matrix{Float64}, k::Int, rng::AbstractRNG)

Initialize both NMF factors with scaled absolute Gaussian values

Arguments

• X::Matrix{Float64}: 2d array of nonnegative floats; the n×p matrix to factorize
• k::Int: The requested factorization rank
• rng::AbstractRNG: Random-number generator used for every draw

Value

A tuple (w, ht) containing an n×k score matrix and a p×k transposed-component matrix. Entries are abs(randn()) * sqrt(mean(X)/k), matching scikit-learn's random NMF initialization while storing Hᵀ for column-major coordinate updates

_nmf_reconstruction_error(X, w, h)

Compute the NMF Frobenius reconstruction error without forming W*H

Arguments

• X: 2d array of floats; the n×p target matrix
• w: 2d array of floats; the n×k score matrix
• h: 2d array of floats; the k×p component matrix

Value

Float; ‖X − W*H‖_F. Accumulates each reconstructed entry directly so it avoids an n×p residual allocation without suffering the cancellation error of the expanded Gram-matrix identity when the reconstruction is nearly exact

_pmd_binary_search(argu, sumabs::Float64)

Find the soft-threshold λ that makes a soft-thresholded vector meet a target L1 budget

Arguments

• argu: 1d array of floats; the vector to be thresholded
• sumabs::Float64: The target ℓ₁/ℓ₂ ratio (the L1 budget)

Value

Float; the threshold λ ≥ 0 such that ‖S(argu,λ)‖₁/‖S(argu,λ)‖₂ ≈ sumabs. Returns 0 when the vector already satisfies the budget. Implements the bisection step of Witten, Tibshirani & Hastie (2009), bracketing λ in [0, max|argu|] and halving for up to 150 iterations or until the bracket is below 1e-6

_pmd_check_v(x::Matrix{Float64}, K::Int)

Initialize the right factors from the leading right singular vectors of the data

Arguments

• x::Matrix{Float64}: 2d array of floats; the (already preprocessed) data matrix
• K::Int: The number of components to initialize

Value

2d array of floats; the p×K matrix of leading right singular vectors of x. Two branches keep the eigenproblem small: for wide data (p > n) it decomposes the n×n Gram matrix x·xᵀ and back-projects, otherwise it decomposes the p×p Gram matrix xᵀ·x directly. SVD initialization is the deterministic starting point recommended in Witten, Tibshirani & Hastie (2009)

_pmd_l1_norm(a)

Compute the ℓ₁/ℓ₂ ratio of a vector

Arguments

• a: 1d array of floats

Value

Float; ‖a‖₁ / ‖a‖₂. This ratio is the quantity the L1 budget constrains; it ranges from 1 (maximally sparse) to √length(a) (fully dense)

_pmd_l1_norm_soft(a, λ::Float64)

Compute the ℓ₁/ℓ₂ ratio of a vector AFTER soft-thresholding it by λ, without materializing the thresholded vector

Arguments

• a: 1d array of floats; the unthresholded vector
• λ::Float64: The soft-threshold to apply

Value

Float; ‖S(a,λ)‖₁ / ‖S(a,λ)‖₂ where S is the soft-threshold operator. Used inside the binary search to test whether a candidate λ meets the L1 budget

_pmd_l1diff(a, b)

Compute the ℓ₁ distance between two vectors

Arguments

• a: 1d array of floats
• b: 1d array of floats

Value

Float; ‖a − b‖₁. Used as the convergence criterion between successive v iterates

_pmd_l2n(a)

Compute the ℓ₂ norm of a vector with a zero-guard

Arguments

• a: 1d array of floats

Value

Float; the Euclidean norm ‖a‖₂, or 0.05 when a is all zeros. The guard mirrors the R package PMA's convention and prevents division by zero in the normalization steps

_pmd_smd!(x, v0, sumabsu::Float64, sumabsv::Float64, niter::Int,
		  u, v, vold, argu, argv)

Compute one rank-1 sparse factor pair by alternating, soft-thresholded power iteration

Arguments

• x: 2d array of floats; the (deflated) data matrix for this component
• v0: 1d array of floats; the initial right factor (from _pmd_check_v)
• sumabsu::Float64: The L1 budget on u (in [1, √n])
• sumabsv::Float64: The L1 budget on v (in [1, √p])
• niter::Int: The maximum number of power-iteration steps
• u: 1d array of floats; preallocated buffer for the left factor, written in place
• v: 1d array of floats; preallocated buffer for the right factor, written in place
• vold: 1d array of floats; preallocated scratch holding the previous v iterate
• argu: 1d array of floats; preallocated scratch of length n
• argv: 1d array of floats; preallocated scratch of length p

Value

Float; the weight d = uᵀ x v for the fitted component. On return u and v hold the unit-ℓ₂, L1-penalized factors. Each step alternates u ← normalize(S(x v, λᵤ)) and v ← normalize(S(xᵀ u, λᵥ)), with λ chosen by _pmd_binary_search to meet the budgets; this is the single-factor PMD(L1,L1) update of Witten, Tibshirani & Hastie (2009). The randn primer on vold only forces the first iteration to run and does not affect the converged fixed point

_pmd_soft(a::Float64, λ::Float64)

Soft-threshold a scalar: shrink toward zero by λ, clamping at zero

Arguments

• a::Float64: The value to threshold
• λ::Float64: The nonnegative threshold

Value

Float; sign(a)·max(|a| − λ, 0). The elementwise solution to the L1-penalized update in Witten, Tibshirani & Hastie (2009)

_pmd_soft_normalize!(out, arg, λ::Float64)

Soft-threshold a vector by λ and rescale it to unit ℓ₂ norm, writing into a preallocated buffer

Arguments

• out: 1d array of floats; the destination buffer, overwritten in place
• arg: 1d array of floats; the vector to threshold
• λ::Float64: The soft-threshold

Value

The buffer out, holding S(arg,λ) normalized to unit ℓ₂ norm (or all zeros, guarded against division by zero). This is the factor update of the PMD power iteration

_prop_var_explained(Xc, V)

Compute the cumulative proportion of variance explained by sparse components

Arguments

• Xc: 2d array of floats; the centered data matrix
• V: 2d array of floats; the p×K matrix of (possibly non-orthonormal) sparse loadings

Value

1d array of floats of length K; entry k is the proportion of total variance explained by the first k loadings jointly. Because sparse loadings need not be orthonormal, the adjusted variance uses the projection ‖Xc·Vₖ(VₖᵀVₖ)⁻¹Vₖᵀ‖²F / ‖Xc‖²F, evaluated here in the equivalent trace form tr((VₖᵀVₖ)⁻¹ VₖᵀXcᵀXcVₖ) / ‖Xc‖²_F. This is the adjusted-variance measure of Shen & Huang / Witten et al. for non-orthogonal sparse PCs

_qnormalize!(P, C)

Normalize the columns of a direction matrix to unit variance under a covariance metric, in place

Arguments

• P: 2d array of floats; the directions (columns), overwritten in place
• C: 2d array of floats; the covariance defining the metric

Value

The matrix P, with each column pⱼ rescaled so that pⱼᵀ C pⱼ = 1 (unit variance of the corresponding canonical variate). Used to normalize the recovered side's directions in the covariance-based CCA solver

_safe_svd!(A)

In-place variant of _safe_svd: compute an SVD, overwriting A, with a robust fallback on convergence failure

Arguments

• A: 2d array of floats; the matrix to decompose, OVERWRITTEN in place by the default path

Value

An SVD factorization object. Tries the in-place svd! first; on a LAPACKException retries with the (non-mutating) QR-iteration SVD. Used in the inner JIVE iterations where the input is scratch that can be destroyed

_safe_svd(A)

Compute an SVD, falling back to a slower but more robust algorithm if the default one fails to converge

Arguments

• A: 2d array of floats; the matrix to decompose

Value

An SVD factorization object. Tries the default divide-and-conquer SVD first; if it raises a LAPACKException (a convergence failure on certain ill-conditioned inputs), retries with the QR-iteration algorithm, which is slower but more stable. Any other error is rethrown

_safe_svdvals(A)

Compute singular values, falling back to a more robust algorithm if the default one fails to converge

Arguments

• A: 2d array of floats; the matrix whose singular values are wanted

Value

A vector of singular values. Like _safe_svd, retries with QR iteration on a LAPACKException and rethrows anything else

_sign_consistency_opt!(V)

Fix the sign of each column of a loading matrix so the largest-magnitude entry is positive, in place

Arguments

• V: 2d array of floats; the loading/direction matrix whose columns are sign-fixed in place

Value

The matrix V, with each column multiplied by ±1 so that its entry of largest absolute value is positive. SVD- and eigen-derived directions are only determined up to a sign, so this canonicalizes the choice to make results reproducible across runs and comparable across implementations. An all-zero column (whose largest entry has sign 0) is left untouched

_soft_threshold_L1!(out, x, nx::Int, absx, ord, ranks)

Soft-threshold a loading vector to keep only its largest-magnitude entries, writing the result into a preallocated buffer

Arguments

• out: 1d array of floats; the destination buffer, overwritten in place
• x: 1d array of floats; the loading vector to threshold
• nx::Int: The number of SMALLEST-magnitude entries to drop (= p − keepX). If nx ≤ 0 the vector is copied through unchanged
• absx: 1d array of floats; preallocated scratch for the entry magnitudes
• ord: 1d array of ints; preallocated scratch for the sort permutation
• ranks: 1d array of ints; preallocated scratch for the magnitude ranks

Value

The buffer out, holding the thresholded vector: the keepX = p − nx largest entries survive (shrunk by λ = the largest dropped magnitude), the rest are zeroed. Implements the L1 penalty Pλ(u) = sign(u)(|u| − λ)₊ of Lê Cao et al. (2011), reparameterized as a keep-count rather than a λ (mixOmics `softthresholding_L1`). The rank computation handles magnitude ties so that tied entries are dropped or kept together

_softcca!(out, a, d)

Soft-threshold a vector into a preallocated buffer

Arguments

• out: 1d array of floats; the destination buffer, overwritten in place
• a: 1d array of floats; the vector to threshold
• d: Float; the nonnegative threshold

Value

The buffer out, holding sign(aᵢ)·max(|aᵢ| − d, 0). The elementwise L1-penalized update of Witten, Tibshirani & Hastie (2009)

_sparse_cca_single_opt!(u, v, x, z, v0, penaltyx, penaltyz, niter,
						vold, zv_buf, xu_buf, argu, argv, su, sv)

Compute one rank-1 sparse canonical vector pair by alternating, soft-thresholded power iteration

Arguments

• u: 1d array of floats; the X canonical vector, written in place
• v: 1d array of floats; the Z canonical vector, written in place
• x: 2d array of floats; the (observations × features) X block for this component
• z: 2d array of floats; the (observations × features) Z block for this component
• v0: 1d array of floats; the initial v (from _fast_init_v or an SVD)
• penaltyx: Float; the X penalty in (0,1], scaled internally to the budget penaltyx·√(#X features)
• penaltyz: Float; the Z penalty in (0,1], scaled to penaltyz·√(#Z features)
• niter: Int; the maximum number of power-iteration steps
• vold: 1d array of floats; scratch holding the previous v iterate
• zv_buf, xu_buf: 1d arrays of floats; scratch for z·v and x·u (sized for the deflated block height)
• argu, argv: 1d arrays of floats; scratch for the raw (pre-threshold) canonical vectors
• su, sv: 1d arrays of floats; scratch for the soft-thresholded vectors

Value

Float; the weight d = uᵀ xᵀ z v = ⟨x u, z v⟩. On return u and v hold the unit-ℓ₂, L1-penalized canonical vectors. Each step alternates u ← normalize(S(xᵀz v, λᵤ)) and v ← normalize(S(zᵀx u, λᵥ)), with the thresholds chosen by _binary_search_opt to meet the penalty budgets — the penalized-CCA update of Witten, Tibshirani & Hastie (2009). The randn primer on vold only forces the first iteration to run and does not affect the converged fixed point

_spca_component!(v, X, c, u, Xv, Xtu, s, vold; tol = 1e-7, niter = 20)

Compute one sparse principal component by soft-thresholded power iteration (deflation variant)

Arguments

• v: 1d array of floats; the right factor, initialized by the caller and overwritten in place with the converged loading
• X: 2d array of floats; the (deflated) data matrix for this component
• c: Float; the L1 budget on v (in [1, √p])
• u: 1d array of floats; preallocated buffer for the left factor (scores direction)
• Xv: 1d array of floats; preallocated scratch of length n
• Xtu: 1d array of floats; preallocated scratch of length p
• s: 1d array of floats; preallocated scratch for _finding_v!
• vold: 1d array of floats; preallocated scratch holding the previous v iterate
• tol: Float; the ℓ₁ convergence tolerance on v. Defaults to 1e-7
• niter::Int: The maximum number of power-iteration steps. Defaults to 20

Value

Float; the component weight d = uᵀ X v. On return v holds the unit-ℓ₂, L1-penalized loading and u the corresponding unit-ℓ₂ scores direction. Each step sets u ← normalize(X v) then v ← project(Xᵀ u) onto the unit ball and L1 budget via _finding_v!, the SPC update of Witten, Tibshirani & Hastie (2009)

_spca_component_orth!(v, X, c, U_prev, u, uold, Xv, Xtu, s, vold, proj;
					 tol = 1e-6, niter = 20)

Compute one sparse principal component whose scores are orthogonal to all previously extracted components

Arguments

• v: 1d array of floats; the right factor, initialized by the caller and overwritten in place with the converged loading
• X: 2d array of floats; the (uncentered-residual) data matrix — note the orthogonal variant decomposes the original Xc, not a deflated residual
• c: Float; the L1 budget on v (in [1, √p])
• U_prev: 2d array of floats; the n×(j-1) matrix of previously extracted scores directions that the new u must be orthogonal to (empty for the first component)
• u: 1d array of floats; preallocated buffer for the left factor, written in place
• uold: 1d array of floats; preallocated scratch holding the previous u iterate
• Xv: 1d array of floats; preallocated scratch of length n
• Xtu: 1d array of floats; preallocated scratch of length p
• s: 1d array of floats; preallocated scratch for _finding_v!
• vold: 1d array of floats; preallocated scratch holding the previous v iterate
• proj: 1d array of floats; preallocated scratch for the orthogonal projection coefficients
• tol: Float; the ℓ₁ convergence tolerance on u and v. Defaults to 1e-6
• niter::Int: The maximum number of power-iteration steps. Defaults to 20

Value

Float; the component weight d = uᵀ X v. Differs from _spca_component! by projecting u onto the orthogonal complement of U_prev each iteration (u ← u − Uprev Uprevᵀ u before normalizing), which enforces orthogonal scores across components instead of relying on deflation — the orthogonal SPC variant of Witten, Tibshirani & Hastie (2009)

_sqdiff(a, b)

Compute the squared ℓ₂ distance between two vectors

Arguments

• a: 1d array of floats
• b: 1d array of floats

Value

Float; Σ(aᵢ − bᵢ)². Used as the convergence criterion between successive loading iterates

_unmap(y::AbstractVector; levels = nothing)

One-hot encode a vector of class labels into a dummy indicator matrix

Arguments

• y::AbstractVector: The class labels, one per observation
• levels: An optional ordering of the classes; when given, it fixes the column order of the output and must list each class exactly once. Defaults to nothing (classes sorted by sort(unique(y)))

Value

A tuple (Yd, classes) where Yd is the n×K one-hot matrix (Yd[i,k] = 1 when observation i is in class k, else 0) and classes is the class labels in column order. This is mixOmics' unmap, the dummy encoding that turns classification into a PLS regression onto the indicator matrix

_validate_nmf_data(X::Matrix{Float64})

Validate a data matrix for nonnegative matrix factorization

Arguments

• X::Matrix{Float64}: 2d array of floats; observations in rows and features in columns

Value

The input matrix X. Throws an ArgumentError if it is empty, contains a non-finite value, or contains a negative value