"""
MoDeVa model wrappers for DirectRS.
MoDirectRSRegressor — post-processes a pre-trained tree ensemble for regression.
MoDirectRSClassifier — binary classification via logistic ridge (IRLS in logit space).
All interpretability methods return ``ValidationResult`` objects following
the MoDeVa standard. Use ``result.plot()`` to visualize and ``result.value``
to access the raw numerical data.
"""
import numpy as np
import pandas as pd
from scipy.special import expit as _sigmoid
from sklearn.base import RegressorMixin, ClassifierMixin
from ..base import ModelBaseRegressor, ModelBaseClassifier
from ...utils.results import ValidationResult
from .core import DirectRSCore
from . import visualize as dviz
class _DirectRSInterpretMixin:
"""Shared interpretability methods for DirectRS models.
All methods call the underlying ``DirectRSCore`` for computation,
build mocharts visualizations via ``visualize.py``, and wrap
everything in ``ValidationResult``.
"""
def get_global_stretch_analysis(self, feature_names=None):
"""Return eigendecomposition of S' matrix.
Returns
-------
ValidationResult
Accessible plots via ``result.plot('eigenvalue_spectrum')``,
``result.plot('feature_activity')``,
``result.plot('feature_loadings')``,
``result.plot('g_matrix')``, or ``result.plot()`` for all.
"""
raw = self._core.get_global_stretch_analysis(feature_names)
names = raw['feature_names']
table = pd.DataFrame({
"Feature": names,
"S'_jj": raw['diag'],
})
options = dviz.viz_stretch_analysis(raw)
return ValidationResult(
key="directrs_stretch_analysis",
model=self.name,
inputs={"feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def explain_local(self, X, feature_names=None, sample_index=0):
"""Exact per-sample additive decomposition in original feature space.
Parameters
----------
X : array-like of shape (n_samples, n_features)
feature_names : list of str, optional
sample_index : int, default=0
Sample to show in the default plot.
Returns
-------
ValidationResult
``result.value`` contains all samples (intercept, contributions
arrays). ``result.plot()`` shows the selected sample.
"""
raw = self._core.explain_local(X, feature_names)
names = raw['feature_names']
c = raw['contributions'][sample_index]
table = pd.DataFrame({
"Name": names,
"Effect": c,
})
options = dviz.viz_local_explanation(raw, sample_index)
return ValidationResult(
key="directrs_local",
model=self.name,
inputs={"sample_index": sample_index,
"feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def importance_global(self, X=None, feature_names=None, mode="slope"):
"""Exact global feature importance in original feature space.
Parameters
----------
X : array-like, optional
Required for data-weighted modes ("contrib_abs", "contrib_rms").
feature_names : list of str, optional
mode : str, default="slope"
"slope", "slope2", "contrib_abs", or "contrib_rms".
Returns
-------
ValidationResult
Horizontal bar plot of feature importance.
"""
raw = self._core.importance_global(X, feature_names, mode)
names = raw['feature_names']
table = pd.DataFrame({
"Name": names,
"Importance": raw['importance'],
})
options = dviz.viz_importance_global(raw)
return ValidationResult(
key="directrs_importance",
model=self.name,
inputs={"mode": mode, "feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def importance_main_interaction(self, X, feature_names=None, n_bins=20):
"""Exact empirical main-effect and interaction decomposition (Eq 10-14).
Parameters
----------
X : array-like of shape (n_samples, n_features)
Data used to compute binned main effects and interaction residuals.
feature_names : list of str, optional
n_bins : int, default=20
Number of quantile bins for continuous features.
Returns
-------
ValidationResult
Grouped bar plot of main vs excess importance per feature.
``result.value`` contains all raw arrays including variance
identity diagnostics (var_f, E_g2, rho, eta2_main, etc.).
"""
raw = self._core.importance_main_interaction(X, feature_names, n_bins)
names = raw['feature_names']
table = pd.DataFrame({
"Feature": names,
"I_main": raw['main'],
"I_excess": raw['excess'],
})
options = dviz.viz_main_interaction(raw)
return ValidationResult(
key="directrs_main_interaction",
model=self.name,
inputs={"n_bins": n_bins, "feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def geometric_interaction_traces(self, feature_names=None, K=3, gamma=0.5):
"""Geometric higher-order interaction tracing via stretch matrix (Section 6).
Parameters
----------
feature_names : list of str, optional
K : int, default=3
Maximum coupling order.
gamma : float, default=0.5
Decay factor for cumulative map (0 < gamma < 1).
Returns
-------
ValidationResult
Accessible plots via ``result.plot('spectrum')``,
``result.plot('adjacency')``,
``result.plot('cumulative_coupling')``,
or ``result.plot()`` for all.
"""
raw = self._core.geometric_interaction_traces(feature_names, K, gamma)
table = pd.DataFrame({
"k": list(range(1, len(raw['T']) + 1)),
"T_k": raw['T'],
"E_k": raw['E'],
})
options = dviz.viz_geometric_traces(raw)
return ValidationResult(
key="directrs_interaction_traces",
model=self.name,
inputs={"K": K, "gamma": gamma,
"feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def get_off_diagonal_analysis(self, feature_names=None, top_k=10):
"""Return top off-diagonal G entries (feature interactions).
Parameters
----------
feature_names : list of str, optional
top_k : int, default=10
Number of top feature pairs to return, ranked by |G_ij|.
Returns
-------
ValidationResult
Horizontal bar plot of top feature-pair coupling strengths.
"""
raw = self._core.get_off_diagonal_analysis(feature_names, top_k)
table = pd.DataFrame(
[(a, b, v) for a, b, v in raw],
columns=["Feature_i", "Feature_j", "Coupling"],
)
options = dviz.viz_off_diagonal(raw)
return ValidationResult(
key="directrs_interactions",
model=self.name,
inputs={"top_k": top_k, "feature_names": feature_names},
value=raw,
table=table,
options=options,
)
def get_feature_importance(self, feature_names=None):
"""Return Ridge-coefficient weighted feature importance.
.. deprecated::
Use ``importance_global()`` instead.
"""
return self._core.get_feature_importance(feature_names)
[docs]
class MoDirectRSRegressor(_DirectRSInterpretMixin, RegressorMixin,
ModelBaseRegressor):
"""DirectRS post-processor for tree ensemble regressors.
Takes a pre-trained tree ensemble (XGBoost, LightGBM, or CatBoost) and builds
a piecewise-linear model that can match or improve the original with per-tree
Ridge regression on geometric embeddings \u03c6_t(x) = [1, S'x].
Parameters
----------
base_model : fitted tree ensemble
Pre-trained XGBoost, LightGBM, or CatBoost regressor (or MoDeVa wrapper).
name : str, optional
Model identifier. Default: "DirectRS".
construction : str, default="C"
Operator construction for geometry extraction:
"A" (gain-weighted), "B" (co-occurrence), "C" (value-weighted), "M" (count).
ridge_alpha : float, default=100.0
Ridge regularization for per-tree linear heads.
n_passes : int, default=1
Number of backfitting coordinate descent passes.
n_trees_used : int, optional
Number of trees to use. Default: all trees.
"""
def __init__(self, base_model, name=None, construction="C",
ridge_alpha=100.0, n_passes=1, n_trees_used=None):
self.base_model = base_model
self.name = name or "DirectRS"
self.construction = construction
self.ridge_alpha = ridge_alpha
self.n_passes = n_passes
self.n_trees_used = n_trees_used
self._core = None
[docs]
def fit(self, X, y, sample_weight=None, X_val=None, y_val=None, verbose=False):
"""Fit DirectRS on training data.
Parameters
----------
X : np.ndarray of shape (n_samples, n_features)
y : np.ndarray of shape (n_samples,)
sample_weight : ignored (kept for API compatibility)
X_val : np.ndarray, optional
y_val : np.ndarray, optional
verbose : bool, default=False
Returns
-------
self
"""
X, y, _ = self._validate_fit_inputs(X, y, sample_weight)
self._core = DirectRSCore(
self.base_model,
construction=self.construction,
ridge_alpha=self.ridge_alpha,
n_passes=self.n_passes,
n_trees_used=self.n_trees_used,
)
self._core.fit(X, y, X_val=X_val, y_val=y_val, verbose=verbose)
return self
def _predict(self, X):
return self._core.predict(X)
[docs]
class MoDirectRSClassifier(_DirectRSInterpretMixin, ClassifierMixin,
ModelBaseClassifier):
"""DirectRS post-processor for tree ensemble classifiers.
Takes a pre-trained tree ensemble classifier and builds a piecewise-linear
model using logistic ridge via IRLS (weighted Ridge on working responses
in logit space). Predictions use sigmoid for probabilities.
Parameters
----------
base_model : fitted tree ensemble
Pre-trained XGBoost, LightGBM, or CatBoost classifier (or MoDeVa wrapper).
name : str, optional
Model identifier. Default: "DirectRS-Cls".
construction : str, default="C"
Operator construction for geometry extraction.
ridge_alpha : float, default=100.0
Ridge regularization for per-tree linear heads.
n_passes : int, default=1
Number of backfitting coordinate descent passes.
n_trees_used : int, optional
Number of trees to use. Default: all trees.
"""
def __init__(self, base_model, name=None, construction="C",
ridge_alpha=100.0, n_passes=1, n_trees_used=None):
self.base_model = base_model
self.name = name or "DirectRS-Cls"
self.construction = construction
self.ridge_alpha = ridge_alpha
self.n_passes = n_passes
self.n_trees_used = n_trees_used
self._core = None
[docs]
def fit(self, X, y, sample_weight=None, X_val=None, y_val=None, verbose=False):
"""Fit DirectRS on training data for classification.
Parameters
----------
X : np.ndarray of shape (n_samples, n_features)
y : np.ndarray of shape (n_samples,)
sample_weight : ignored
X_val : np.ndarray, optional
y_val : np.ndarray, optional
verbose : bool, default=False
Returns
-------
self
"""
X, y, _ = self._validate_fit_inputs(X, y, sample_weight)
self._core = DirectRSCore(
self.base_model,
construction=self.construction,
ridge_alpha=self.ridge_alpha,
n_passes=self.n_passes,
n_trees_used=self.n_trees_used,
)
self._core.fit(X, y, X_val=X_val, y_val=y_val, verbose=verbose,
classification=True)
return self
def _predict_proba(self, X):
"""Predict probabilities via sigmoid on logit-space output."""
raw = self._core.predict(X)
prob_1 = _sigmoid(raw)
return np.column_stack([1.0 - prob_1, prob_1])