pod.compute_kriging_loo_residuals

compute_kriging_loo_residuals(gpr, X_2d, y)

Computes exact LOO predictions, variances, standardized residuals e_i, and outlier scaling factor gamma according to Malkiel et al. (2026).

Inverts the augmented covariance matrix S = [[K + alpha*I, F], [F^T, 0]] to obtain B = S^-1, calculating Leave-One-Out means mu_{-i} and variances sigma_{-i}^2. It then derives standardized residuals e_i = (y_i - mu_{-i}) / sigma_{-i} and outlier scale factor gamma = max(1.0, max|e_i| / 3.0).

Parameters

Name	Type	Description	Default
gpr	GaussianProcessRegressor	A fitted scikit-learn Gaussian Process model.	required
X_2d	np.ndarray	2D matrix of input training coordinates (N_samples, N_features).	required
y	np.ndarray	1D array of observed responses (N_samples,).	required

Returns

Name	Type	Description
	Tuple[np.ndarray, np.ndarray, np.ndarray, float]	Tuple[np.ndarray, np.ndarray, np.ndarray, float]: - loo_means: Array of Leave-One-Out predicted mean responses. - loo_stds: Array of Leave-One-Out predicted standard deviations. - std_residuals: Array of standardized LOO residuals e_i = (y_i - mu_{-i}) / sigma_{-i}. - gamma: Outlier scaling factor gamma = max(1.0, max|e_i| / 3.0).

Examples

import numpy as np
from sklearn.gaussian_process import GaussianProcessRegressor
from digiqual.pod import compute_kriging_loo_residuals

X = np.linspace(0, 5, 20).reshape(-1, 1)
y = 2.0 * X.flatten() + np.random.normal(0, 0.1, 20)

gpr = GaussianProcessRegressor()
gpr.fit(X, y)

loo_means, loo_stds, std_res, gamma = compute_kriging_loo_residuals(gpr, X, y)
print(f"Outlier scaling factor gamma: {gamma:.3f}")