Principal Covariates Regression (PCovR)

PCovR

class skmatter.decomposition.PCovR(mixing=0.5, n_components=None, svd_solver='auto', tol=1e-12, space='auto', regressor=None, iterated_power='auto', random_state=None)

Bases: _BasePCA, LinearModel

Principal Covariates Regression, as described in [deJong1992] determines a latent-space projection \(\mathbf{T}\) which minimizes a combined loss in supervised and unsupervised tasks.

This projection is determined by the eigendecomposition of a modified gram matrix \(\mathbf{\tilde{K}}\)

\[\mathbf{\tilde{K}} = \alpha \mathbf{X} \mathbf{X}^T + (1 - \alpha) \mathbf{\hat{Y}}\mathbf{\hat{Y}}^T\]

where \(\alpha\) is a mixing parameter and \(\mathbf{X}\) and \(\mathbf{\hat{Y}}\) are matrices of shapes \((n_{samples}, n_{features})\) and \((n_{samples}, n_{properties})\), respectively, which contain the input and approximate targets. For \((n_{samples} < n_{features})\), this can be more efficiently computed using the eigendecomposition of a modified covariance matrix \(\mathbf{\tilde{C}}\)

\[\mathbf{\tilde{C}} = \alpha \mathbf{X}^T \mathbf{X} + (1 - \alpha) \left(\left(\mathbf{X}^T \mathbf{X}\right)^{-\frac{1}{2}} \mathbf{X}^T \mathbf{\hat{Y}}\mathbf{\hat{Y}}^T \mathbf{X} \left(\mathbf{X}^T \mathbf{X}\right)^{-\frac{1}{2}}\right)\]

For all PCovR methods, it is strongly suggested that \(\mathbf{X}\) and \(\mathbf{Y}\) are centered and scaled to unit variance, otherwise the results will change drastically near \(\alpha \to 0\) and \(\alpha \to 0\). This can be done with the companion preprocessing classes, where

>>> from skmatter.preprocessing import StandardFlexibleScaler as SFS
>>>
>>> # Set column_wise to True when the columns are relative to one another,
>>> # False otherwise.
>>> scaler = SFS(column_wise=True)
>>>
>>> scaler.fit(A) # replace with your matrix
>>> A = scaler.transform(A)

Parameters:

• mixing (float, default=0.5) – mixing parameter, as described in PCovR as \({\alpha}\), here named to avoid confusion with regularization parameter alpha

• n_components (int, float or str, default=None) –

  Number of components to keep. if n_components is not set all components are kept:

  n_components == min(n_samples, n_features)
  

• svd_solver ({'auto', 'full', 'arpack', 'randomized'}, default='auto') –

  If auto :

  The solver is selected by a default policy based on X.shape and n_components: if the input data is larger than 500x500 and the number of components to extract is lower than 80% of the smallest dimension of the data, then the more efficient ‘randomized’ method is enabled. Otherwise the exact full SVD is computed and optionally truncated afterwards.

  If full :

  run exact full SVD calling the standard LAPACK solver via scipy.linalg.svd and select the components by postprocessing

  If arpack :

  run SVD truncated to n_components calling ARPACK solver via scipy.sparse.linalg.svds. It requires strictly 0 < n_components < min(X.shape)

  If randomized :

  run randomized SVD by the method of Halko et al.

• tol (float, default=1e-12) – Tolerance for singular values computed by svd_solver == ‘arpack’. Must be of range [0.0, infinity).

• space ({'feature', 'sample', 'auto'}, default='auto') – whether to compute the PCovR in sample or feature space default=`sample` when \({n_{samples} < n_{features}}\) and feature when \({n_{features} < n_{samples}}\)

• regressor ({Ridge, RidgeCV, LinearRegression, precomputed}, default=None) – regressor for computing approximated \({\mathbf{\hat{Y}}}\). The regressor should be one sklearn.linear_model.Ridge, sklearn.linear_model.RidgeCV, or sklearn.linear_model.LinearRegression. If a pre-fitted regressor is provided, it is used to compute \({\mathbf{\hat{Y}}}\). Note that any pre-fitting of the regressor will be lost if PCovR is within a composite estimator that enforces cloning, e.g., sklearn.compose.TransformedTargetRegressor or sklearn.pipeline.Pipeline with model caching. In such cases, the regressor will be re-fitted on the same training data as the composite estimator. If precomputed, we assume that the y passed to the fit function is the regressed form of the targets \({\mathbf{\hat{Y}}}\). If None, sklearn.linear_model.Ridge(‘alpha’:1e-6, ‘fit_intercept’:False, ‘tol’:1e-12) is used as the regressor.

• iterated_power (int or 'auto', default='auto') – Number of iterations for the power method computed by svd_solver == ‘randomized’. Must be of range [0, infinity).

• random_state (int, RandomState instance or None, default=None) – Used when the ‘arpack’ or ‘randomized’ solvers are used. Pass an int for reproducible results across multiple function calls.

mixing

mixing parameter, as described in PCovR as \({\alpha}\)

Type:

float, default=0.5

tol

Tolerance for singular values computed by svd_solver == ‘arpack’. Must be of range [0.0, infinity).

Type:

float, default=1e-12

space

whether to compute the PCovR in sample or feature space default=`sample` when \({n_{samples} < n_{features}}\) and feature when \({n_{features} < n_{samples}}\)

Type:

{‘feature’, ‘sample’, ‘auto’}, default=’auto’

n_components

The estimated number of components, which equals the parameter n_components, or the lesser value of n_features and n_samples if n_components is None.

Type:

int

pxt_

the projector, or weights, from the input space \(\mathbf{X}\) to the latent-space projection \(\mathbf{T}\)

Type:

ndarray of size \(({n_{samples}, n_{components}})\)

pty_

the projector, or weights, from the latent-space projection \(\mathbf{T}\) to the properties \(\mathbf{Y}\)

Type:

ndarray of size \(({n_{components}, n_{properties}})\)

pxy_

the projector, or weights, from the input space \(\mathbf{X}\) to the properties \(\mathbf{Y}\)

Type:

ndarray of size \(({n_{samples}, n_{properties}})\)

explained_variance_

The amount of variance explained by each of the selected components.

Equal to n_components largest eigenvalues of the PCovR-modified covariance matrix of \(\mathbf{X}\).

Type:

ndarray of shape (n_components,)

singular_values_

The singular values corresponding to each of the selected components.

Type:

ndarray of shape (n_components,)

Examples

>>> import numpy as np
>>> from skmatter.decomposition import PCovR
>>> X = np.array([[-1, 1, -3, 1], [1, -2, 1, 2], [-2, 0, -2, -2], [1, 0, 2, -1]])
>>> Y = np.array([[ 0, -5], [-1, 1], [1, -5], [-3, 2]])
>>> pcovr = PCovR(mixing=0.1, n_components=2)
>>> pcovr.fit(X, Y)
PCovR(mixing=0.1, n_components=2, space='sample')
>>> T = pcovr.transform(X)
    [[ 3.2630561 ,  0.06663787],
     [-2.69395511, -0.41582771],
     [ 3.48683147, -0.83164387],
     [-4.05593245,  1.18083371]]
>>> Yp = pcovr.predict(X)
    [[ 0.01374656, -5.00943466],
     [-1.02804032,  1.06737777],
     [ 0.98167556, -4.9830631 ],
     [-2.99627428,  1.98241962]]

fit(X, Y, W=None)

Fit the model with X and Y. Depending on the dimensions of X, calls either _fit_feature_space or _fit_sample_space

Parameters:

• X (ndarray, shape (n_samples, n_features)) –

  Training data, where n_samples is the number of samples and n_features is the number of features.

  It is suggested that \(\mathbf{X}\) be centered by its column- means and scaled. If features are related, the matrix should be scaled to have unit variance, otherwise \(\mathbf{X}\) should be scaled so that each feature has a variance of 1 / n_features.

• Y (ndarray, shape (n_samples, n_properties)) –

  Training data, where n_samples is the number of samples and n_properties is the number of properties

  It is suggested that \(\mathbf{X}\) be centered by its column- means and scaled. If features are related, the matrix should be scaled to have unit variance, otherwise \(\mathbf{Y}\) should be scaled so that each feature has a variance of 1 / n_features.

  If the passed regressor = precomputed, it is assumed that Y is the regressed form of the properties, \({\mathbf{\hat{Y}}}\).

• W (ndarray, shape (n_features, n_properties)) – Regression weights, optional when regressor=`precomputed`. If not passed, it is assumed that W = np.linalg.lstsq(X, Y, self.tol)[0]

_fit_feature_space(X, Y, Yhat)

In feature-space PCovR, the projectors are determined by:

\[\mathbf{\tilde{C}} = \alpha \mathbf{X}^T \mathbf{X} + (1 - \alpha) \left(\left(\mathbf{X}^T \mathbf{X}\right)^{-\frac{1}{2}} \mathbf{X}^T \mathbf{\hat{Y}}\mathbf{\hat{Y}}^T \mathbf{X} \left(\mathbf{X}^T \mathbf{X}\right)^{-\frac{1}{2}}\right)\]

where

\[\mathbf{P}_{XT} = (\mathbf{X}^T \mathbf{X})^{-\frac{1}{2}} \mathbf{U}_\mathbf{\tilde{C}}^T \mathbf{\Lambda}_\mathbf{\tilde{C}}^{\frac{1}{2}}\]

\[\mathbf{P}_{TX} = \mathbf{\Lambda}_\mathbf{\tilde{C}}^{-\frac{1}{2}} \mathbf{U}_\mathbf{\tilde{C}}^T (\mathbf{X}^T \mathbf{X})^{\frac{1}{2}}\]

\[\mathbf{P}_{TY} = \mathbf{\Lambda}_\mathbf{\tilde{C}}^{-\frac{1}{2}} \mathbf{U}_\mathbf{\tilde{C}}^T (\mathbf{X}^T \mathbf{X})^{-\frac{1}{2}} \mathbf{X}^T \mathbf{Y}\]

_fit_sample_space(X, Y, Yhat, W)

In sample-space PCovR, the projectors are determined by:

\[\mathbf{\tilde{K}} = \alpha \mathbf{X} \mathbf{X}^T + (1 - \alpha) \mathbf{\hat{Y}}\mathbf{\hat{Y}}^T\]

where

\[\mathbf{P}_{XT} = \left(\alpha \mathbf{X}^T + (1 - \alpha) \mathbf{W} \mathbf{\hat{Y}}^T\right) \mathbf{U}_\mathbf{\tilde{K}} \mathbf{\Lambda}_\mathbf{\tilde{K}}^{-\frac{1}{2}}\]

\[\mathbf{P}_{TX} = \mathbf{\Lambda}_\mathbf{\tilde{K}}^{-\frac{1}{2}} \mathbf{U}_\mathbf{\tilde{K}}^T \mathbf{X}\]

\[\mathbf{P}_{TY} = \mathbf{\Lambda}_\mathbf{\tilde{K}}^{-\frac{1}{2}} \mathbf{U}_\mathbf{\tilde{K}}^T \mathbf{Y}\]

transform(X=None)

Apply dimensionality reduction to X.

X is projected on the first principal components as determined by the modified PCovR distances.

Parameters:

X (ndarray, shape (n_samples, n_features)) – New data, where n_samples is the number of samples and n_features is the number of features.

predict(X=None, T=None)

Predicts the property values using regression on X or T

inverse_transform(T)

Transform data back to its original space.

\[\mathbf{\hat{X}} = \mathbf{T} \mathbf{P}_{TX} = \mathbf{X} \mathbf{P}_{XT} \mathbf{P}_{TX}\]

Parameters:

T (ndarray, shape (n_samples, n_components)) – Projected data, where n_samples is the number of samples and n_components is the number of components.

Returns:

X_original ndarray, shape (n_samples, n_features)

score(X, Y, T=None)

Return the (negative) total reconstruction error for X and Y, defined as:

\[\ell_{X} = \frac{\lVert \mathbf{X} - \mathbf{T}\mathbf{P}_{TX} \rVert ^ 2}{\lVert \mathbf{X}\rVert ^ 2}\]

and

\[\ell_{Y} = \frac{\lVert \mathbf{Y} - \mathbf{T}\mathbf{P}_{TY} \rVert ^ 2}{\lVert \mathbf{Y}\rVert ^ 2}\]

The negative loss \(-\ell = -(\ell_{X} + \ell{Y})\) is returned for easier use in sklearn pipelines, e.g., a grid search, where methods named ‘score’ are meant to be maximized.

Parameters:

• X (ndarray of shape (n_samples, n_features)) – The data.

• Y (ndarray of shape (n_samples, n_properties)) – The target.

Returns:

loss (float) – Negative sum of the loss in reconstructing X from the latent-space projection T and the loss in predicting Y from the latent-space projection T

Kernel PCovR

class skmatter.decomposition.KernelPCovR(mixing=0.5, n_components=None, svd_solver='auto', regressor=None, kernel='linear', gamma=None, degree=3, coef0=1, kernel_params=None, center=False, fit_inverse_transform=False, tol=1e-12, n_jobs=None, iterated_power='auto', random_state=None)

Bases: _BasePCA, LinearModel

Kernel Principal Covariates Regression, as described in [Helfrecht2020] determines a latent-space projection \(\mathbf{T}\) which minimizes a combined loss in supervised and unsupervised tasks in the reproducing kernel Hilbert space (RKHS).

This projection is determined by the eigendecomposition of a modified gram matrix \(\mathbf{\tilde{K}}\)

\[\mathbf{\tilde{K}} = \alpha \mathbf{K} + (1 - \alpha) \mathbf{\hat{Y}}\mathbf{\hat{Y}}^T\]

where \(\alpha\) is a mixing parameter, \(\mathbf{K}\) is the input kernel of shape \((n_{samples}, n_{samples})\) and \(\mathbf{\hat{Y}}\) is the target matrix of shape \((n_{samples}, n_{properties})\).

Parameters:

• mixing (float, default=0.5) – mixing parameter, as described in PCovR as \({\\alpha}\)

• n_components (int, float or str, default=None) –

  Number of components to keep. if n_components is not set all components are kept:

  n_components == n_samples
  

• svd_solver ({'auto', 'full', 'arpack', 'randomized'}, default='auto') –

  If auto :

  The solver is selected by a default policy based on X.shape and n_components: if the input data is larger than 500x500 and the number of components to extract is lower than 80% of the smallest dimension of the data, then the more efficient ‘randomized’ method is enabled. Otherwise the exact full SVD is computed and optionally truncated afterwards.

  If full :

  run exact full SVD calling the standard LAPACK solver via scipy.linalg.svd and select the components by postprocessing

  If arpack :

  run SVD truncated to n_components calling ARPACK solver via scipy.sparse.linalg.svds. It requires strictly 0 < n_components < min(X.shape)

  If randomized :

  run randomized SVD by the method of Halko et al.

• regressor ({instance of sklearn.kernel_ridge.KernelRidge, precomputed, None}, default=None) –

  The regressor to use for computing the property predictions \(\\hat{\\mathbf{Y}}\). A pre-fitted regressor may be provided. If the regressor is not None, its kernel parameters (kernel, gamma, degree, coef0, and kernel_params) must be identical to those passed directly to KernelPCovR.

  If precomputed, we assume that the y passed to the fit function is the regressed form of the targets \({\mathbf{\hat{Y}}}\).

kernel: “linear” | “poly” | “rbf” | “sigmoid” | “cosine” | “precomputed”

Kernel. Default=”linear”.

gamma: float, default=None

Kernel coefficient for rbf, poly and sigmoid kernels. Ignored by other kernels.

degree: int, default=3

Degree for poly kernels. Ignored by other kernels.

coef0: float, default=1

Independent term in poly and sigmoid kernels. Ignored by other kernels.

kernel_params: mapping of str to any, default=None

Parameters (keyword arguments) and values for kernel passed as callable object. Ignored by other kernels.

center: bool, default=False

Whether to center any computed kernels

fit_inverse_transform: bool, default=False

Learn the inverse transform for non-precomputed kernels. (i.e. learn to find the pre-image of a point)

tol: float, default=1e-12

Tolerance for singular values computed by svd_solver == ‘arpack’ and for matrix inversions. Must be of range [0.0, infinity).

n_jobs: int, default=None

The number of parallel jobs to run. None means 1 unless in a joblib.parallel_backend context. -1 means using all processors.

iterated_powerint or ‘auto’, default=’auto’

Number of iterations for the power method computed by svd_solver == ‘randomized’. Must be of range [0, infinity).

random_stateint, RandomState instance or None, default=None

Used when the ‘arpack’ or ‘randomized’ solvers are used. Pass an int for reproducible results across multiple function calls.

pt__

pseudo-inverse of the latent-space projection, which can be used to contruct projectors from latent-space

Type:

ndarray of size \(({n_{components}, n_{components}})\)

pkt_

the projector, or weights, from the input kernel \(\\mathbf{K}\) to the latent-space projection \(\\mathbf{T}\)

Type:

ndarray of size \(({n_{samples}, n_{components}})\)

pky_

the projector, or weights, from the input kernel \(\\mathbf{K}\) to the properties \(\\mathbf{Y}\)

Type:

ndarray of size \(({n_{samples}, n_{properties}})\)

pty_

the projector, or weights, from the latent-space projection \(\\mathbf{T}\) to the properties \(\\mathbf{Y}\)

Type:

ndarray of size \(({n_{components}, n_{properties}})\)

ptx_

the projector, or weights, from the latent-space projection \(\\mathbf{T}\) to the feature matrix \(\\mathbf{X}\)

Type:

ndarray of size \(({n_{components}, n_{features}})\)

X_fit_

The data used to fit the model. This attribute is used to build kernels from new data.

Type:

ndarray of shape (n_samples, n_features)

Examples

>>> import numpy as np
>>> from skmatter.decomposition import KernelPCovR
>>> from skmatter.preprocessing import StandardFlexibleScaler as SFS
>>> from sklearn.kernel_ridge import KernelRidge
>>>
>>> X = np.array([[-1, 1, -3, 1], [1, -2, 1, 2], [-2, 0, -2, -2], [1, 0, 2, -1]])
>>> X = SFS().fit_transform(X)
>>> Y = np.array([[ 0, -5], [-1, 1], [1, -5], [-3, 2]])
>>> Y = SFS(column_wise=True).fit_transform(Y)
>>>
>>> kpcovr = KernelPCovR(mixing=0.1, n_components=2, regressor=KernelRidge(kernel='rbf', gamma=1), kernel='rbf', gamma=1)
>>> kpcovr.fit(X, Y)
    KernelPCovR(gamma=1, kernel='rbf', mixing=0.1, n_components=2,
                regressor=KernelRidge(gamma=1, kernel='rbf'))
>>> T = kpcovr.transform(X)
    [[-0.61261285, -0.18937908],
     [ 0.45242098,  0.25453465],
     [-0.77871824,  0.04847559],
     [ 0.91186937, -0.21211816]]
>>> Yp = kpcovr.predict(X)
    [[ 0.5100212 , -0.99488463],
     [-0.18992219,  0.82064368],
     [ 1.11923584, -1.04798016],
     [-1.5635827 ,  1.11078662]]
>>> kpcovr.score(X, Y)
    -0.520388347837897

fit(X, Y, W=None)

Fit the model with X and Y.

Parameters:

• X (ndarray, shape (n_samples, n_features)) –

  Training data, where n_samples is the number of samples and n_features is the number of features.

  It is suggested that \(\mathbf{X}\) be centered by its column- means and scaled. If features are related, the matrix should be scaled to have unit variance, otherwise \(\mathbf{X}\) should be scaled so that each feature has a variance of 1 / n_features.

• Y (ndarray, shape (n_samples, n_properties)) –

  Training data, where n_samples is the number of samples and n_properties is the number of properties

  It is suggested that \(\mathbf{X}\) be centered by its column- means and scaled. If features are related, the matrix should be scaled to have unit variance, otherwise \(\mathbf{Y}\) should be scaled so that each feature has a variance of 1 / n_features.

• W (ndarray, shape (n_samples, n_properties)) – Regression weights, optional when regressor=`precomputed`. If not passed, it is assumed that W = np.linalg.lstsq(K, Y, self.tol)[0]

Returns:

self (object) – Returns the instance itself.

transform(X)

Apply dimensionality reduction to X.

X is projected on the first principal components as determined by the modified Kernel PCovR distances.

Parameters:

X (ndarray, shape (n_samples, n_features)) – New data, where n_samples is the number of samples and n_features is the number of features.

predict(X=None)

Predicts the property values

inverse_transform(T)

Transform input data back to its original space.

\[\mathbf{\hat{X}} = \mathbf{T} \mathbf{P}_{TX} = \mathbf{K} \mathbf{P}_{KT} \mathbf{P}_{TX}\]

Similar to KPCA, the original features are not always recoverable, as the projection is computed from the kernel features, not the original features, and the mapping between the original and kernel features is not one-to-one.

Parameters:

T (ndarray, shape (n_samples, n_components)) – Projected data, where n_samples is the number of samples and n_components is the number of components.

Returns:

X_original ndarray, shape (n_samples, n_features)

score(X, Y)

Computes the (negative) loss values for KernelPCovR on the given predictor and response variables. The loss in \(\mathbf{K}\), as explained in [Helfrecht2020] does not correspond to a traditional Gram loss \(\mathbf{K} - \mathbf{TT}^T\). Indicating the kernel between set A and B as \(\mathbf{K}_{AB}\), the projection of set A as \(\mathbf{T}_A\), and with N and V as the train and validation/test set, one obtains

\[\ell=\frac{\operatorname{Tr}\left[\mathbf{K}_{VV} - 2 \mathbf{K}_{VN} \mathbf{T}_N (\mathbf{T}_N^T \mathbf{T}_N)^{-1} \mathbf{T}_V^T +\mathbf{T}_V(\mathbf{T}_N^T \mathbf{T}_N)^{-1} \mathbf{T}_N^T \mathbf{K}_{NN} \mathbf{T}_N (\mathbf{T}_N^T \mathbf{T}_N)^{-1} \mathbf{T}_V^T\right]}{\operatorname{Tr}(\mathbf{K}_{VV})}\]

The negative loss is returned for easier use in sklearn pipelines, e.g., a grid search, where methods named ‘score’ are meant to be maximized.

Parameters:

• X (independent (predictor) variable) –

• Y (dependent (response) variable) –

Returns:

L (Negative sum of the KPCA and KRR losses, with the KPCA loss) – determined by the reconstruction of the kernel