"""MoDeVa wrappers for the fuseKernel fused kernel-ridge model.
``MoFuseKernelRegressor`` and ``MoFuseKernelClassifier`` make fuseKernel a first-class MoDeVa
model. They train on a ``DataSet``, plug into ``ModelZoo`` and the whole ``TestSuite``, and add:
* **Inherent FANOVA interpretation** -- ``ts.interpret_global_fi/ei/effect`` and
``ts.interpret_local_fi/ei`` work through the standard contract (``interpret``,
``modeva_effects_``, ``predict_effect``). The main effects and pairwise interactions are
fuseKernel's functional-ANOVA decomposition of the fused predictor; interaction candidates are
ranked by the spectral channel's ARD relevance. Needs ``use_spectral=True``.
* **Inherent GP predictive intervals** (regressor) -- ``predict_interval`` returns fuseKernel's
closed-form Gaussian-process posterior interval instead of conformal, so ``diagnose_reliability``
and ``model.predict_interval`` reflect the model's own uncertainty.
* **Channel decomposition** -- ``channel_contributions`` returns the exact additive per-channel
split of each prediction as a ``ValidationResult``.
All interpretability methods return ``ValidationResult`` objects following the MoDeVa standard.
"""
import numpy as np
import pandas as pd
from sklearn.base import RegressorMixin, ClassifierMixin
from ..base import ModelBaseRegressor, ModelBaseClassifier
from ...utils.results import ValidationResult
from ._impl import FuseKernelModel, FuseConfig, SpectralInterpreter
from ._impl.interpret import channel_contributions as _fk_channel_contributions
# ---------------------------------------------------------------- effect adapter
class _FusedEffectAdapter:
"""Presents the fused fuseKernel model to a :class:`SpectralInterpreter`: ARD relevance comes
from the spectral channel's MS-SKM (``smix_``), but predictions are the *fused* model's, so the
partial-dependence / functional-ANOVA effects decompose the whole fused predictor."""
def __init__(self, fk_model):
self._fk = fk_model
smix = None
for ch in getattr(fk_model, "channels_", []):
if getattr(ch, "name", None) == "spectral":
smix = ch.sk.m.smix_
break
if smix is None:
raise ValueError("inherent FANOVA interpretation needs a spectral channel; "
"construct the model with use_spectral=True.")
self.smix_ = smix
self.is_clf_ = fk_model.is_clf_
self.feature_names_in_ = None
self.scaler_ = None
def predict(self, X):
return self._fk.predict(np.asarray(X, dtype=np.float64))
def predict_proba(self, X):
return self._fk.predict_proba(np.asarray(X, dtype=np.float64))
# ---------------------------------------------------------------- FANOVA + extras mixin
class _FuseKernelInterpretMixin:
"""Inherent FANOVA interpretation (the MoDeVa ``interpret`` contract) plus fuseKernel-specific
interpretability returned as ``ValidationResult`` objects. Shared by regressor and classifier."""
# -- the MoDeVa FANOVA contract -------------------------------------------------
def interpret(self, dataset):
"""Build the inherent FANOVA interpreter. Precomputes fuseKernel's main-effect curves and
pairwise-interaction surfaces over a reference sample, then exposes them through the
standard :class:`InterpretFANOVA` so ``ts.interpret_*`` work. Needs a spectral channel."""
from ...testsuite.interpret.fanova.base import InterpretFANOVA
names = list(dataset.feature_names)
d = len(names)
X_ref = np.asarray(dataset.get_X_y_data(dataset="train")[0], dtype=np.float64)
if len(X_ref) > self.interpret_ref_size:
rng = np.random.RandomState(self.random_state)
X_ref = X_ref[np.sort(rng.choice(len(X_ref), self.interpret_ref_size, replace=False))]
si = SpectralInterpreter(_FusedEffectAdapter(self.estimator_), feature_names=names)
self._interp_si_ = si
# main effects: one centered partial-dependence curve per feature
self._eff_main_ = {}
for j in range(d):
gj, ej = si.main_effect(j, X_ref, n_grid=self.interpret_n_grid)
self._eff_main_[j] = (gj, ej)
# interactions: pairs among the most ARD-relevant features, capped
imp = si.feature_importance(normalize=False)
pool = [int(j) for j in np.argsort(imp)[::-1][:self.max_interaction_features]]
pairs = [(pool[a], pool[b]) for a in range(len(pool)) for b in range(a + 1, len(pool))]
pairs = pairs[:self.max_interactions]
ngi = max(8, self.interpret_n_grid // 2)
self._eff_inter_ = {}
interaction = {}
for idx, (j, k) in enumerate(pairs):
gj, gk, surf = si.interaction_effect(j, k, X_ref, n_grid=ngi)
self._eff_inter_[idx] = (j, k, gj, gk, surf)
interaction[f"{names[j]} & {names[k]}"] = {"fidx": (j, k), "interaction_idx": idx}
self.modeva_effects_ = {
"main_effect": {names[j]: {"fidx": (j,)} for j in range(d)},
"interaction": interaction,
}
self.modeva_intercept_ = float(si._response_batched(X_ref, 4096).mean())
return InterpretFANOVA(model=self, dataset=dataset)
def predict_main_effect(self, X):
"""Per-feature main-effect raw predictions, shape (n, n_features)."""
X = np.asarray(X, dtype=np.float64)
out = np.zeros((X.shape[0], len(self._eff_main_)))
for j, (gj, ej) in self._eff_main_.items():
out[:, j] = np.interp(X[:, j], gj, ej)
return out
def predict_interaction(self, X):
"""Pairwise-interaction raw predictions, shape (n, n_interactions)."""
X = np.asarray(X, dtype=np.float64)
out = np.zeros((X.shape[0], len(self._eff_inter_)))
for idx, (j, k, gj, gk, surf) in self._eff_inter_.items():
out[:, idx] = SpectralInterpreter._bilinear(gj, gk, surf, X[:, j], X[:, k])
return out
def predict_effect(self, fidx, X):
"""Raw prediction of one main effect (len-1 fidx) or pairwise interaction (len-2 fidx)."""
if isinstance(fidx, (int, np.integer)):
fidx = (int(fidx),)
fidx = tuple(int(i) for i in fidx)
if len(fidx) == 1:
return self.predict_main_effect(X)[:, fidx[0]].ravel()
for _, item in self.modeva_effects_.get("interaction", {}).items():
if tuple(item["fidx"]) == fidx:
return self.predict_interaction(X)[:, item["interaction_idx"]].ravel()
return np.zeros(np.asarray(X).shape[0])
# -- fuseKernel-specific: exact per-channel decomposition -----------------------
def channel_contributions(self, X):
"""Exact additive per-channel decomposition of the fused prediction (regression).
Returns a ``ValidationResult``: ``value`` holds the per-channel contribution arrays and the
intercept; ``table`` is the mean absolute contribution per channel; ``plot()`` shows the bar.
"""
import mocharts as mc
contribs, intercept = _fk_channel_contributions(self.estimator_, np.asarray(X, dtype=np.float64))
names = list(contribs.keys())
mean_abs = [float(np.mean(np.abs(v))) for v in contribs.values()]
table = (pd.DataFrame({"Channel": names, "MeanAbsContribution": mean_abs})
.sort_values("MeanAbsContribution", ascending=False).reset_index(drop=True))
options = mc.barplot(x=names, y=mean_abs, orient="horizontal")
options.set_title("fuseKernel channel contributions")
options.set_xaxis(axis_name="Mean |contribution|")
options.set_tooltip(precision=4)
return ValidationResult(
key="fusekernel_channel_contributions",
model=self.name,
inputs=None,
value={"contributions": contribs, "intercept": intercept},
table=table,
options=options.render(),
)
# -- fuseKernel-specific: weak-cluster diagnosis and repair ---------------------
def _require_fitted(self):
if getattr(self, "estimator_", None) is None:
raise RuntimeError("fit the model before calling this method.")
def diagnose_weak_clusters(self, dataset, n_clusters: int = 5, **kw):
"""Per-cluster train/test performance breakdown of the fitted fused kernel.
Nyström spectral clustering of the model's own kernel partitions the data into
``n_clusters`` regions; the model metric is reported per region on train and
test, so regions where the model underperforms (large train/test gap or low
headline metric) are exposed. Works for both the classic two-channel and the
general (spectral / multi-depth) paths.
Returns a ``ValidationResult``: ``table`` is the per-cluster metric breakdown
(with an ``ALL`` aggregate row); ``value`` holds the cluster labels, spectral
embeddings and the weakest-cluster ranking; ``plot()`` shows the per-cluster
test metric as a bar.
"""
import mocharts as mc
from ._impl.diagnostics import _primary
self._require_fitted()
X_tr, y_tr, _ = dataset.get_X_y_data(dataset="train")
X_te, y_te, _ = dataset.get_X_y_data(dataset="test")
X_tr = np.asarray(X_tr, dtype=np.float64); y_tr = np.asarray(y_tr).ravel()
X_te = np.asarray(X_te, dtype=np.float64); y_te = np.asarray(y_te).ravel()
res = self.estimator_.diagnose(X_tr, y_tr, X_te, y_te, n_clusters=n_clusters, **kw)
prim, higher_better = _primary(res.task)
per = res.table[res.table["cluster"] != "ALL"]
metric_col = f"test_{prim}" if f"test_{prim}" in per.columns else f"train_{prim}"
clusters = [str(c) for c in per["cluster"].tolist()]
scores = [float(v) for v in per[metric_col].tolist()]
options = mc.barplot(x=clusters, y=scores, orient="vertical")
options.set_title(f"FuseKernel per-cluster {metric_col}")
options.set_xaxis(axis_name="cluster")
options.set_yaxis(axis_name=prim)
options.set_tooltip(precision=4)
return ValidationResult(
key="fusekernel_diagnose_weak_clusters",
model=self.name,
inputs={"n_clusters": n_clusters, "metric": prim},
value={
"labels_train": res.labels_train,
"labels_test": res.labels_test,
"embedding_train": res.embedding_train,
"embedding_test": res.embedding_test,
"worst_clusters": res.worst_clusters(k=min(3, n_clusters)),
},
table=res.table,
options=options.render(),
)
#: valid MS-SKM decode solvers for the spectral channel
_SOLVERS = ("auto", "lanczos", "matfree", "nystrom")
def _build_config(self, task):
# surface the spectral-channel decode solver as a first-class option; an explicit
# key in ``spectral_params`` still wins, otherwise fall back to ``self.solver``.
if self.solver is not None and self.solver not in self._SOLVERS:
raise ValueError(
f"solver must be one of {self._SOLVERS}, got {self.solver!r}")
spectral_params = dict(self.spectral_params) if self.spectral_params else {}
if self.solver is not None:
spectral_params.setdefault("solver", self.solver)
return FuseConfig(
task=task, backend=self.backend, gbdt_params=self.gbdt_params,
use_xgb=self.use_xgb, use_rbf=self.use_rbf, use_spectral=self.use_spectral,
tree_depths=self.tree_depths, spectral_params=spectral_params or None,
fit_method=self.fit_method, residual_nw=self.residual_nw, seed=self.random_state,
)
# ============================================================================ regressor
[docs]
class MoFuseKernelRegressor(_FuseKernelInterpretMixin, RegressorMixin, ModelBaseRegressor):
"""fuseKernel fused kernel-ridge regressor for MoDeVa.
Parameters
----------
name : str, optional
use_xgb, use_rbf, use_spectral : bool
Which kernel channels to fuse (defaults: tree + spectral). ``use_xgb`` turns on the
tree co-membership channel regardless of ``backend``.
tree_depths : tuple of int, optional
One co-membership kernel per depth (multi-depth tree fusion).
spectral_params : dict, optional
MS-SKM kwargs (H / K / kernel / ...). None -> the fuseKernel defaults
(``kernel="laplace", H=4, K=8, solver="nystrom"``).
fit_method : {"grid", "adam", "nlml", "oof", "gcv", "sure"}, default="grid"
Fusion-weight selection (grid/adam are leakage-free, query-scored).
solver : {"nystrom", "auto", "lanczos", "matfree"}, default="nystrom"
How the spectral channel decodes its kernel (only used when
``use_spectral=True``). ``"nystrom"`` is the linear-in-n low-rank decode
and the fastest/most scalable default; ``"lanczos"`` is the exact dense
decode (best for small data); ``"auto"`` uses dense below ~20k rows and
switches to Nystrom above; ``"matfree"`` is a matrix-free CG solve. An
explicit ``"solver"`` key in ``spectral_params`` overrides this.
backend : {"xgboost", "lightgbm", "catboost"}, default="xgboost"
Gradient-boosted ensemble that defines the leaf co-membership partition.
gbdt_params : dict, optional
Native params for the chosen backend. None -> that backend's defaults.
residual_nw : passthrough fuseKernel option.
max_interaction_features, max_interactions : interaction screening for ``interpret``.
interpret_n_grid, interpret_ref_size : fANOVA grid / reference-sample size.
random_state : int
"""
def __init__(self, name=None, use_xgb=True, use_rbf=False, use_spectral=True,
tree_depths=None, spectral_params=None, fit_method="grid",
solver="nystrom", backend="xgboost", gbdt_params=None, residual_nw=False,
max_interaction_features=6, max_interactions=10,
interpret_n_grid=30, interpret_ref_size=1000, random_state=0):
self.name = name or "FuseKernel"
self.use_xgb = use_xgb
self.use_rbf = use_rbf
self.use_spectral = use_spectral
self.tree_depths = tree_depths
self.spectral_params = spectral_params
self.fit_method = fit_method
self.solver = solver
self.backend = backend
self.gbdt_params = gbdt_params
self.residual_nw = residual_nw
self.max_interaction_features = max_interaction_features
self.max_interactions = max_interactions
self.interpret_n_grid = interpret_n_grid
self.interpret_ref_size = interpret_ref_size
self.random_state = random_state
self.estimator_ = None
self.feature_names_ = None
def fit(self, X, y, sample_weight=None, feature_names=None):
X, y, _ = self._validate_fit_inputs(X, y, sample_weight)
self.feature_names_ = (list(feature_names) if feature_names is not None
else [f"x{j}" for j in range(X.shape[1])])
self.estimator_ = FuseKernelModel(self._build_config("regression")).fit(X, y)
return self
def _predict(self, X):
return self.estimator_.predict(np.asarray(X, dtype=np.float64))
# -- inherent GP predictive intervals (override the conformal default) ----------
[docs]
def calibrate_interval(self, X, y, alpha=0.1, max_depth=5):
"""fuseKernel's intervals are the GP posterior -- no fitting needed; just record alpha."""
self._gp_alpha_ = float(alpha)
self.calibration_qval_ = 0.0 # mark "calibrated"; predict_interval is overridden below
[docs]
def predict_interval(self, X):
"""Closed-form GP posterior prediction interval at level ``1 - alpha``."""
from scipy.stats import norm
alpha = getattr(self, "_gp_alpha_", 0.1)
mean, var = self.estimator_.predict_dist(np.asarray(X, dtype=np.float64))
z = float(norm.ppf(1 - alpha / 2))
half = z * np.sqrt(var)
return np.vstack([mean - half, mean + half]).T
[docs]
def predict_dist(self, X):
"""GP predictive ``(mean, variance)`` in target units."""
return self.estimator_.predict_dist(np.asarray(X, dtype=np.float64))
# ============================================================================ classifier
[docs]
class MoFuseKernelClassifier(_FuseKernelInterpretMixin, ClassifierMixin, ModelBaseClassifier):
"""fuseKernel fused kernel-ridge classifier for MoDeVa (binary and multiclass).
The fused KRR decodes one-hot targets to per-class scores; ``predict_proba`` is a
temperature-calibrated softmax. Same parameters as :class:`MoFuseKernelRegressor`.
"""
def __init__(self, name=None, use_xgb=True, use_rbf=False, use_spectral=True,
tree_depths=None, spectral_params=None, fit_method="grid",
solver="nystrom", backend="xgboost", gbdt_params=None,
max_interaction_features=6, max_interactions=10,
interpret_n_grid=30, interpret_ref_size=1000, random_state=0):
self.name = name or "FuseKernel-Cls"
self.use_xgb = use_xgb
self.use_rbf = use_rbf
self.use_spectral = use_spectral
self.tree_depths = tree_depths
self.spectral_params = spectral_params
self.fit_method = fit_method
self.solver = solver
self.backend = backend
self.gbdt_params = gbdt_params
self.residual_nw = False
self.max_interaction_features = max_interaction_features
self.max_interactions = max_interactions
self.interpret_n_grid = interpret_n_grid
self.interpret_ref_size = interpret_ref_size
self.random_state = random_state
self.estimator_ = None
self.feature_names_ = None
def fit(self, X, y, sample_weight=None, feature_names=None):
y_raw = np.asarray(y)
X, _, _ = self._validate_fit_inputs(X, y, sample_weight) # sets classes_/label_binarizer_
self.feature_names_ = (list(feature_names) if feature_names is not None
else [f"x{j}" for j in range(X.shape[1])])
self.estimator_ = FuseKernelModel(self._build_config("classification")).fit(X, y_raw)
return self
[docs]
def predict_proba(self, X):
return self.estimator_.predict_proba(np.asarray(X, dtype=np.float64))
def _predict_proba(self, X):
return self.predict_proba(X)
[docs]
def predict(self, X):
return self.estimator_.predict(np.asarray(X, dtype=np.float64))
def _predict(self, X):
return self.predict(X)
[docs]
def decision_function(self, X):
proba = self.predict_proba(X)
return proba[:, 1] if proba.shape[1] == 2 else proba