Skip to content

GP Methods

The cca_zoo.gp module provides GaussianProcessCCA, a nonlinear multiview CCA method that uses a Gaussian process with a joint (non-additive) kernel as the per-view encoder. It has no optional dependency: the kernel machinery is built entirely from scikit-learn's own GaussianProcessRegressor, RBF, ConstantKernel and KernelCenterer, with scipy.linalg doing the fit, all already required by cca_zoo.


Background

GaussianProcessCCA minimises the same unconstrained Eckart-Young (EY) objective used by the stochastic *_EY models in cca_zoo.linear, by DCCAEY in cca_zoo.deep, by TreeCCA in cca_zoo.tree, and by GAMCCA in cca_zoo.gam (all share the exact same implementation, in cca_zoo._utils._ey):

\[ \mathcal{L}_{EY} = -2 \operatorname{tr}(C) + \operatorname{tr}(V V) \]

where, for embeddings \(Z_i = f_i(X_i)\), \(C\) is the mean pairwise cross-covariance (including \(i = j\) terms) and \(V\) the mean auto-covariance across all views. Where GAMCCA sums one univariate B-spline term per input feature — additive, and so structurally unable to represent an interaction between two features of the same view — GaussianProcessCCA fits a Gaussian process with a joint RBF kernel (one length scale shared by all features, following the data's spread by default) directly over each view's raw, joint feature vector: a genuine, non-additive function of all of that view's features at once. As a Gaussian process it also comes with a posterior standard deviation for every score for free.

Each encoder writes \(f_i(x) = k_c(x, U_i)^\top B_i\) for a fixed kernel \(k\), centred in feature space, and a fixed set of basis ("inducing") points \(U_i\) — every training row by default, or n_inducing of them selected via sklearn.cluster.kmeans_plusplus for larger datasets. The coefficients minimise the EY loss plus the RKHS-norm penalty \(\tfrac12 \alpha_i \operatorname{tr}(B_i^\top K_i B_i)\). Writing \(B_i\) in the eigenbasis of the inducing points' kernel turns the encoder into a linear map of their Nyström features and the penalty into a ridge, so the global minimiser over every view jointly is a generalized eigenproblem, solved in closed form: no iterations, no initialisation, no local optima. With a linear kernel and no penalty it is exactly CCA. Kernel hyperparameters are fixed — pass kernel= explicitly, or tune it externally (e.g. with cca_zoo.model_selection.GridSearchCV).

When to use: Nonlinear multiview CCA where cross-view structure depends on a genuine interaction between two or more features of the same view (e.g. \(x_1 x_2\)), which an additive GAM cannot represent and which a tree ensemble can only approximate via multivariate splits (and, at default settings, may be starved of joint feature access — see TreeCCA's colsample_bytree). A GP with a joint kernel represents such an interaction directly. When the true relationship is smooth and additive, GAMCCA is typically a cheaper and equally accurate choice; a GP over the raw feature vector scales cubically in the number of training samples (exact inference), so expect fitting to be markedly slower than GAMCCA or TreeCCA on large datasets unless n_inducing is set.


Basic usage

from cca_zoo.gp import GaussianProcessCCA

model = GaussianProcessCCA(n_components=1).fit([X1, X2])
z1, z2 = model.transform([X1, X2])
corr = model.score([X1, X2])  # mean canonical correlation

# GaussianProcessCCA also supports more than two views
model3 = GaussianProcessCCA(n_components=1).fit([X1, X2, X3])

Predictive uncertainty

Because each per-component encoder is a Gaussian process, posterior_std gives each latent component's posterior standard deviation, the uncertainty of transform's scores:

means, stds = model.transform([X1, X2]), model.posterior_std([X1, X2])

stds[i] has the same shape as means[i] ((n_samples, n_components)) and is the posterior standard deviation of view i's latent component — larger away from the training data, smaller near it, exactly as for any other GP posterior. This does not depend on the fitted coefficients at all (a standard GP fact: posterior variance only involves the kernel, the noise level, and the design points), so it is computed under the (uncentred) GP prior implied by the same kernel, noise level, and inducing points as the mean fit.

Read it as epistemic uncertainty about where the encoder is determined by the training data, not as a measure of local noise or prediction difficulty: it tracks training density, so where measurement noise is largest relative to the signal in the dense region, the most confident predictions can be the least accurate. Its magnitude is on the kernel's prior scale, not the encoder output's scale, so use it to rank or compare points rather than as a calibrated interval.

GaussianProcessCCA has no linear weight matrices and no per-feature decomposition analogous to GAMCCA's shape_function (the kernel is not additive across features). feature_importances_per_view_ is therefore permutation-based: the mean squared change in each view's latent scores when a feature's training values are shuffled, normalised to sum to 1 per view.


Scaling to larger datasets: the sparse (inducing-point) approximation

Exact inference costs \(O(n^3)\) per fit, using every training row as a basis point — impractical somewhere in the low thousands of samples. Setting n_inducing switches each encoder to a reduced-rank ("subset of regressors") construction (Quiñonero-Candela & Rasmussen, 2005), using only n_inducing basis points — an actual, well-spread subset of the training rows, selected via sklearn.cluster.kmeans_plusplus's seeding — reducing fitting to \(O(n \, m^2 + m^3)\) for \(m\) basis points instead of \(O(n^3)\):

model = GaussianProcessCCA(n_components=1, n_inducing=200).fit(
    [X1, X2]
)  # X1, X2 have many samples

n_inducing values at or above the number of training samples are equivalent to (and internally fall back to) exact inference. Larger n_inducing trades speed for a closer approximation to the exact posterior; there's no universal default, since how many inducing points are "enough" depends on how smooth/low-rank the true underlying function is — start with a few hundred and check whether increasing it changes the held-out canonical correlation.


Key parameters

Parameter Description
kernel Fixed kernel, cloned per view, or one per view. None (default) uses ConstantKernel(1.0) * RBF(sqrt(p * X.var() / 2)) for each view's own feature count p, the length scale of scikit-learn's gamma="scale". Pass an RBF with a per-feature length_scale array for an ARD kernel.
alpha Ridge (RKHS-norm) penalty strength, also used as the noise level for the posterior-variance calculation.
n_inducing Number of basis ("inducing") points (see above). None (default) uses every training row (exact inference).
random_state Seed for selecting inducing points when n_inducing is set.

Practical notes

  • GaussianProcessCCA supports 2 or more views.
  • Exact inference (n_inducing=None) is \(O(n^3)\) in the number of training samples; set n_inducing for datasets beyond a few thousand samples, or consider GAMCCA (if the relationship is additive) or TreeCCA instead.
  • No optional dependency is required (unlike cca_zoo.tree, which needs xgboost/lightgbm): GaussianProcessCCA is built entirely on scikit-learn's GaussianProcessRegressor, RBF, ConstantKernel, KernelCenterer and kmeans_plusplus.