Skip to content

cca_zoo.gam

Spline-based nonlinear CCA methods: generalized additive models with mgcv-style P-spline smooths (GAMCCA) and multivariate adaptive regression splines in the style of earth (MARSCCA). Both fit every view jointly by a closed-form penalised eigenproblem. Built on scikit-learn's SplineTransformer and SciPy — no optional dependency required.


GAMCCA

GAMCCA(
    n_components: int = 1,
    *,
    center: bool = True,
    k: int | list[int] = 10,
    m: int
    | tuple[int, int]
    | list[int | tuple[int, int]] = 2,
    sp: float | list[float] = 0.01,
)

Bases: BaseModel

Nonlinear CCA with generalized additive model encoders.

Each view's encoder is \(f_i(x) = \sum_j s_{ij}(x_j)\), one mgcv P-spline smooth s(x_j, bs="ps", k=k, m=m) per feature, fitted to minimise the EY loss (:mod:cca_zoo._utils._ey) plus the smoothing penalty \(\tfrac12 \mathrm{sp}_i \sum_j \beta_{ij}^\top D^\top D \beta_{ij}\), with \(D\) the m[1]-th order difference matrix. The fit is a single generalized eigenproblem, solved in closed form at its global optimum. mgcv's GCV and REML have no EY counterpart, so choose sp by cross-validation. An additive model cannot represent within-view interactions; see :class:~cca_zoo.gam.MARSCCA.

Parameters:

Name Type Description Default
n_components int

Number of latent dimensions. Default is 1.

1
center bool

Whether to centre each view. Default is True.

True
k int | list[int]

B-splines per feature, as mgcv's k. Per-view. Default is 10.

10
m int | tuple[int, int] | list[int | tuple[int, int]]

(order, penalty_order), as mgcv's m: splines of degree order + 1 and a penalty_order-th difference penalty; an int sets both. Per-view. Default is 2.

2
sp float | list[float]

Smoothing parameter; larger is smoother. Per-view. Default is 0.01.

0.01

Attributes:

Name Type Description
encoders_

Fitted per-view encoders.

References

Chapman, J., Wang, H.-T., Wells, L., & Wiesner, J. (2021). CCA-Zoo: A collection of Regularized, Deep Learning based, Kernel, and Probabilistic CCA methods in a scikit-learn style framework. Journal of Open Source Software, 6(68), 3823. Chapman, J., Wells, L., & Lawry Aguila, A. (2024). Unconstrained Stochastic CCA: Unifying Multiview and Self-Supervised Learning. arXiv:2310.01012. Wood, S. N. (2017). Generalized Additive Models: An Introduction with R (2nd ed.). Chapman and Hall/CRC.

Eilers, P. H., & Marx, B. D. (1996). Flexible smoothing with B-splines and penalties. Statistical Science, 11(2), 89-121.

Examples:

>>> import numpy as np
>>> from cca_zoo.gam import GAMCCA
>>> rng = np.random.default_rng(0)
>>> X1 = rng.standard_normal((200, 5))
>>> X2 = rng.standard_normal((200, 5))
>>> model = GAMCCA(n_components=2, k=[8, 12], sp=[0.1, 10.0]).fit([X1, X2])
>>> model.shape_function(view=0, feature=0, x=X1[:, 0]).shape
(200, 2)
Source code in cca_zoo/gam/_gamcca.py
def __init__(
    self,
    n_components: int = 1,
    *,
    center: bool = True,
    k: int | list[int] = 10,
    m: int | tuple[int, int] | list[int | tuple[int, int]] = 2,
    sp: float | list[float] = 0.01,
) -> None:
    super().__init__(n_components=n_components, center=center)
    self.k = k
    self.m = m
    self.sp = sp

fit

fit(views: list[ArrayLike], y: None = None) -> GAMCCA

Fit the model.

Parameters:

Name Type Description Default
views list[ArrayLike]

Arrays of shape (n_samples, n_features_i), one per view.

required
y None

Ignored.

None

Returns:

Type Description
GAMCCA

self.

Raises:

Type Description
ValueError

If k is too small for a view's spline order.

Source code in cca_zoo/gam/_gamcca.py
def fit(self, views: list[ArrayLike], y: None = None) -> GAMCCA:
    """Fit the model.

    Args:
        views: Arrays of shape (n_samples, n_features_i), one per view.
        y: Ignored.

    Returns:
        self.

    Raises:
        ValueError: If ``k`` is too small for a view's spline order.
    """
    views_ = self._setup_fit(views)
    m_views = self.n_views_
    k_ = perview_parameter("k", self.k, 10, m_views)
    m_: list[int | tuple[int, int]] = perview_parameter("m", self.m, 2, m_views)
    sp_ = perview_parameter("sp", self.sp, 0.01, m_views)
    for i, (k, m) in enumerate(zip(k_, m_)):
        order, penalty_order = _spline_orders(m)
        if k < order + 2 or penalty_order >= k:
            raise ValueError(
                f"k={k} is too small for m={m} in view {i}: a P-spline of "
                f"order {order} needs k >= {order + 2}, and the penalty "
                f"order must be below k."
            )
    encoders = [_GamEncoder(X, k, m) for X, k, m in zip(views_, k_, m_)]
    raw_bases = [enc.raw_basis(X) for enc, X in zip(encoders, views_)]
    # The stacked Gram of the centred, constrained bases, from the sparse
    # raw bases: Z'(S - 1 mu')'(S - 1 mu')Z = Z'(S'S - n mu mu')Z.
    stacked = sparse.hstack(raw_bases).tocsr()
    mean = np.concatenate([enc.basis_mean_ for enc in encoders])
    n = stacked.shape[0]
    scale = m_views * (n - 1)
    raw_gram = _sparse_gram(stacked) - n * np.outer(mean, mean)
    # H G H for the block-diagonal reflection H = I - 2UU' is a rank-2p
    # update of G, O(d^2 p) rather than two dense d^3 products.
    reflectors = scipy.linalg.block_diag(*[enc.reflectors_ for enc in encoders])
    gram_reflectors = raw_gram @ reflectors
    reflected = (
        raw_gram
        - 2 * reflectors @ gram_reflectors.T
        - 2 * gram_reflectors @ reflectors.T
        + 4 * reflectors @ (reflectors.T @ gram_reflectors) @ reflectors.T
    )
    free = np.concatenate([enc.free_ for enc in encoders])
    gram = reflected[np.ix_(free, free)] / scale
    view = np.repeat(
        np.arange(m_views), [enc.constraint_.shape[1] for enc in encoders]
    )
    reduced = [
        full_rank_reparametrisation(
            enc.basis_operator(raw, np.sqrt(scale)),
            gram[np.ix_(view == i, view == i)],
            np.sqrt(sp) * enc.penalty_factor_ @ enc.constraint_,
            np.sqrt(sp) * enc.penalty_norm_,
        )
        for i, (enc, raw, sp) in enumerate(zip(encoders, raw_bases, sp_))
    ]
    rows = [row_i for _, row_i, _ in reduced]
    reduced_gram = np.block(
        [
            [
                rows[i].T @ gram[np.ix_(view == i, view == j)] @ rows[j]
                for j in range(m_views)
            ]
            for i in range(m_views)
        ]
    )
    reduced_view = np.repeat(np.arange(m_views), [r.shape[1] for r in rows])
    coefficients = penalised_gram_ey_closed_form(
        reduced_gram,
        reduced_view,
        self.n_components,
        [penalty for _, _, penalty in reduced],
    )
    for enc, (lift_i, _, _), coef in zip(encoders, reduced, coefficients):
        enc.coef_ = enc.constraint_ @ lift_i @ coef

    self.encoders_: list[_GamEncoder] = encoders
    self._fit_maps_and_importances(views_)
    return self

shape_function

shape_function(
    view: int, feature: int, x: ArrayLike
) -> np.ndarray

Evaluate one feature's fitted smooth, plot.gam's partial effect.

Parameters:

Name Type Description Default
view int

Index of the view.

required
feature int

Index of the feature.

required
x ArrayLike

Raw values of that feature, shape (n,).

required

Returns:

Type Description
ndarray

The feature's contribution, shape (n, n_components).

Source code in cca_zoo/gam/_gamcca.py
def shape_function(self, view: int, feature: int, x: ArrayLike) -> np.ndarray:
    """Evaluate one feature's fitted smooth, ``plot.gam``'s partial effect.

    Args:
        view: Index of the view.
        feature: Index of the feature.
        x: Raw values of that feature, shape (n,).

    Returns:
        The feature's contribution, shape (n, n_components).
    """
    check_is_fitted(self)
    x_arr = np.asarray(x, dtype=float) - self.means_[view][feature]
    return self.encoders_[view].feature_term(feature, x_arr)

MARSCCA

MARSCCA(
    n_components: int = 1,
    *,
    center: bool = True,
    degree: int | list[int] = 1,
    nk: int | list[int | None] | None = None,
    thresh: float = 0.001,
    minspan: int | list[int | None] | None = None,
    endspan: int | list[int | None] | None = None,
    alpha: float | list[float] = 0.01,
    nprune: int | None = None,
)

Bases: BaseModel

Nonlinear CCA with multivariate adaptive regression spline encoders.

Each view's encoder is a MARS model (Friedman, 1991), a linear combination of products of up to degree hinges \(\max(0, \pm(x_j - t))\), fitted to minimise the ridge-penalised EY loss. As in R's earth, a forward pass grows each basis by the hinge pair that best absorbs the EY gradient, refitting every view jointly in closed form, and a backward pass prunes it to nprune terms. earth's GCV has no EY counterpart, so choose nprune by cross-validation.

Parameters:

Name Type Description Default
n_components int

Number of latent dimensions. Default is 1.

1
center bool

Whether to centre each view. Default is True.

True
degree int | list[int]

Maximum hinges per term; 1 is additive, 2 allows pairwise interactions. Per-view. Default is 1.

1
nk int | list[int | None] | None

Maximum terms per view in the forward pass; None is earth's min(200, max(20, 2 * n_features)) less its intercept. Per-view. Default is None.

None
thresh float

Stop the forward pass once a round lowers the loss by less than this fraction. Default is 0.001.

0.001
minspan int | list[int | None] | None

Minimum support points between knots; 0 is Friedman's rule, None widens it to at most 20 knots per feature. Per-view. Default is None.

None
endspan int | list[int | None] | None

Support points at either end that may not carry a knot, doubled for interactions; None is Friedman's rule. Per-view. Default is None.

None
alpha float | list[float]

Ridge penalty on the coefficient of every basis function at unit variance, so that it ignores each feature's units. Per-view. Default is 0.01.

0.01
nprune int | None

Total terms, across views, kept by the backward pass; None keeps them all. Default is None.

None

Attributes:

Name Type Description
encoders_

Fitted per-view encoders; encoders_[i].coef_ has one row per entry of :meth:basis_functions.

backward_path_

(view, term) in the order the backward pass deleted them.

backward_loss_

Training EY loss after each number of deletions.

n_removed_

Terms the backward pass pruned, the first n_removed_ of backward_path_; zero when nprune is None.

References

Chapman, J., Wang, H.-T., Wells, L., & Wiesner, J. (2021). CCA-Zoo: A collection of Regularized, Deep Learning based, Kernel, and Probabilistic CCA methods in a scikit-learn style framework. Journal of Open Source Software, 6(68), 3823. Chapman, J., Wells, L., & Lawry Aguila, A. (2024). Unconstrained Stochastic CCA: Unifying Multiview and Self-Supervised Learning. arXiv:2310.01012. Friedman, J. H. (1991). Multivariate Adaptive Regression Splines. The Annals of Statistics, 19(1), 1-67.

Examples:

>>> import numpy as np
>>> from cca_zoo.gam import MARSCCA
>>> rng = np.random.default_rng(0)
>>> X1 = rng.standard_normal((200, 5))
>>> X2 = rng.standard_normal((200, 5))
>>> model = MARSCCA(degree=2, nk=[10, 20]).fit([X1, X2])
>>> len(model.basis_functions(0)) <= 10
True
Source code in cca_zoo/gam/_marscca.py
def __init__(
    self,
    n_components: int = 1,
    *,
    center: bool = True,
    degree: int | list[int] = 1,
    nk: int | list[int | None] | None = None,
    thresh: float = 0.001,
    minspan: int | list[int | None] | None = None,
    endspan: int | list[int | None] | None = None,
    alpha: float | list[float] = 0.01,
    nprune: int | None = None,
) -> None:
    super().__init__(n_components=n_components, center=center)
    self.degree = degree
    self.nk = nk
    self.thresh = thresh
    self.minspan = minspan
    self.endspan = endspan
    self.alpha = alpha
    self.nprune = nprune

fit

fit(views: list[ArrayLike], y: None = None) -> MARSCCA

Fit the model: earth's forward pass, backward pass and final refit.

Parameters:

Name Type Description Default
views list[ArrayLike]

Arrays of shape (n_samples, n_features_i), one per view.

required
y None

Ignored.

None

Returns:

Type Description
MARSCCA

self.

Raises:

Type Description
ValueError

If nprune is below the number of views, or a view has no admissible knot.

Source code in cca_zoo/gam/_marscca.py
def fit(self, views: list[ArrayLike], y: None = None) -> MARSCCA:
    """Fit the model: earth's forward pass, backward pass and final refit.

    Args:
        views: Arrays of shape (n_samples, n_features_i), one per view.
        y: Ignored.

    Returns:
        self.

    Raises:
        ValueError: If ``nprune`` is below the number of views, or a view
            has no admissible knot.
    """
    views_ = self._setup_fit(views)
    m = self.n_views_
    if self.nprune is not None and self.nprune < m:
        raise ValueError(
            f"nprune={self.nprune} is fewer than the number of views "
            f"({m}); every view keeps at least one term."
        )
    k = self.n_components
    alpha_ = perview_parameter("alpha", self.alpha, 0.01, m)
    terms, raw_bases = self._forward_pass(views_, alpha_)

    removed, path_loss = _backward_path(_standardised(raw_bases), k, alpha_)
    view_of = np.repeat(np.arange(m), [len(t) for t in terms])
    column_of = np.concatenate([np.arange(len(t)) for t in terms])
    self.backward_path_: list[tuple[int, _Term]] = [
        (int(view_of[c]), terms[view_of[c]][column_of[c]]) for c in removed
    ]
    self.backward_loss_: np.ndarray = path_loss
    n_removed = 0 if self.nprune is None else max(len(view_of) - self.nprune, 0)
    gone = set(removed[:n_removed].tolist())
    keep = [
        np.array([c not in gone for c in np.flatnonzero(view_of == i)])
        for i in range(m)
    ]
    terms = [
        [t for t, kept in zip(ts, mask) if kept] for ts, mask in zip(terms, keep)
    ]
    raw_bases = [raw[:, mask] for raw, mask in zip(raw_bases, keep)]
    self.n_removed_: int = n_removed

    coefficients = penalised_basis_ey_closed_form(
        _standardised(raw_bases), k, alpha_
    )
    self.encoders_: list[_MarsEncoder] = [
        _MarsEncoder(t, raw.mean(axis=0), coef / raw.std(axis=0)[:, None])
        for t, raw, coef in zip(terms, raw_bases, coefficients)
    ]
    self._fit_maps_and_importances(views_)
    return self

basis_functions

basis_functions(view: int) -> list[str]

Readable form of one view's selected basis functions.

Knots are in raw feature units, e.g. "h(x3 - 0.52) * h(1.1 - x0)" with h(u) = max(0, u).

Parameters:

Name Type Description Default
view int

Index of the view.

required

Returns:

Type Description
list[str]

One string per basis function, matching the rows of

list[str]

encoders_[view].coef_.

Source code in cca_zoo/gam/_marscca.py
def basis_functions(self, view: int) -> list[str]:
    """Readable form of one view's selected basis functions.

    Knots are in raw feature units, e.g. ``"h(x3 - 0.52) * h(1.1 - x0)"``
    with ``h(u) = max(0, u)``.

    Args:
        view: Index of the view.

    Returns:
        One string per basis function, matching the rows of
        ``encoders_[view].coef_``.
    """
    check_is_fitted(self)
    means = self.means_[view]

    def factor(feature: int, knot: float, sign: int) -> str:
        raw_knot = knot + means[feature]
        if sign < 0:
            return f"h({raw_knot:.4g} - x{feature})"
        op = "-" if raw_knot >= 0 else "+"
        return f"h(x{feature} {op} {abs(raw_knot):.4g})"

    return [
        " * ".join(factor(*f) for f in term) for term in self.encoders_[view].terms_
    ]