From be7fa01cd1b3b39415a8ecdba517c30566a395b6 Mon Sep 17 00:00:00 2001 From: ChrisW09 <50968720+ChrisW09@users.noreply.github.com> Date: Fri, 16 Oct 2026 14:20:01 +0200 Subject: [PATCH] fix(splines): drop absolute clipping thresholds from the M-spline basis MSplineTransformer zeroed every basis function whose knot span was below 1e-6, clipped the values to 1e6 and the whole design to +/-1e3. M-spline values scale like 1 / span, so these absolute limits made the basis depend on the units of the feature: on small-scale features, or with tight knots, the peaks were flattened at 1000 and the integrals fell below one, and on tiny-range features the output became all zero. The clips are removed, so each basis function integrates to one at any feature scale and rescaling x by c rescales the basis by 1 / c. Only a basis function whose support is degenerate relative to the knot range (zero width, possible only for a repeated knot) is returned as zeros. Fixes #63 --- pretab/expansion/spline/m_spline.py | 21 +++++----- .../spline/test_spline_expansions.py | 39 +++++++++++++++++++ 2 files changed, 51 insertions(+), 9 deletions(-) diff --git a/pretab/expansion/spline/m_spline.py b/pretab/expansion/spline/m_spline.py index f98682e..beac100 100644 --- a/pretab/expansion/spline/m_spline.py +++ b/pretab/expansion/spline/m_spline.py @@ -74,7 +74,14 @@ def _feature_suffix(self) -> str: return "ms" def _mspline_basis(self, x: np.ndarray, knots: np.ndarray, basis_idx: int) -> np.ndarray: - """Compute a single M-spline basis function.""" + """Compute a single M-spline basis function. + + ``M_k = (p + 1) / (t_{k+p+1} - t_k) * B_k``, which integrates to one. The + values scale like ``1 / span`` and are left unclipped, so the basis does + not depend on the units of the feature. Only a basis function whose + support is degenerate relative to the knot range (zero width, possible + only for a knot vector with a repeated knot) is returned as zeros. + """ n_coef = len(knots) - self.degree - 1 coef = np.zeros(n_coef) coef[basis_idx] = 1.0 @@ -82,17 +89,13 @@ def _mspline_basis(self, x: np.ndarray, knots: np.ndarray, basis_idx: int) -> np values = np.nan_to_num(spline(x), nan=0.0) knot_span = knots[basis_idx + self.degree + 1] - knots[basis_idx] - if knot_span > 1e-6: - values = values * (self.degree + 1) / knot_span - values = np.clip(values, 0.0, 1e6) - else: - values = np.zeros(len(x)) - return np.maximum(values, 0.0) + if knot_span <= np.finfo(float).eps * (knots[-1] - knots[0]): + return np.zeros(len(x)) + return np.maximum(values * (self.degree + 1) / knot_span, 0.0) def _design_matrix(self, x: np.ndarray, knots: np.ndarray) -> np.ndarray: n_coef = len(knots) - self.degree - 1 design = np.zeros((len(x), n_coef)) for i in range(n_coef): design[:, i] = self._mspline_basis(x, knots, i) - design = np.nan_to_num(design, nan=0.0, posinf=0.0, neginf=0.0) - return np.clip(design, -1e3, 1e3) + return np.nan_to_num(design, nan=0.0) diff --git a/tests/expansion/spline/test_spline_expansions.py b/tests/expansion/spline/test_spline_expansions.py index a55d9ca..347607f 100644 --- a/tests/expansion/spline/test_spline_expansions.py +++ b/tests/expansion/spline/test_spline_expansions.py @@ -116,6 +116,45 @@ def test_mspline_handles_nan(): assert np.isfinite(Xt).all() +@pytest.mark.parametrize("scale", [1.0, 1e-2, 1e-3, 1e-7]) +def test_mspline_basis_integrates_to_one_at_any_feature_scale(scale): + """Regression guard for issue #63: absolute clipping thresholds flattened the + basis on small-scale features and zeroed it on tiny-range ones.""" + from scipy.integrate import trapezoid + + X = np.random.default_rng(0).uniform(0, 1, (2000, 1)) * scale + transformer = MSplineTransformer(output_dim=8, placement_strategy="uniform").fit(X) + grid = np.linspace(X.min(), X.max(), 200_001) + integrals = trapezoid(transformer.transform(grid.reshape(-1, 1)), grid, axis=0) + np.testing.assert_allclose(integrals, 1.0, atol=1e-3) + + +def test_mspline_basis_is_scale_equivariant(): + X = np.random.default_rng(0).lognormal(0, 2, (2000, 1)) + unit = MSplineTransformer(output_dim=8).fit(X).transform(X) + scaled = MSplineTransformer(output_dim=8).fit(X * 1e-4).transform(X * 1e-4) + # M-spline values scale like 1 / span, so rescaling x by c rescales M by 1 / c. + np.testing.assert_allclose(scaled * 1e-4, unit, rtol=1e-8) + + +def test_mspline_basis_is_not_clipped_after_the_preprocessor_scaling(): + """Tight quantile knots on a skewed feature, after the default min-max scaling to + [-1, 1], pushed M-spline peaks above the absolute clip at 1000.""" + import pandas as pd + from scipy.integrate import trapezoid + + from pretab import Preprocessor + + x = np.random.default_rng(0).lognormal(0, 2, 2000) + pre = Preprocessor(numerical_method="mspline", target_aware=False, placement_strategy="quantile") + pre.fit(pd.DataFrame({"x": x})) + grid = np.linspace(x.min(), x.max(), 200_001) + basis = np.asarray(pre.transform(pd.DataFrame({"x": grid}))) + scaled = (grid - x.min()) / (x.max() - x.min()) * 2.0 - 1.0 + assert basis.max() > 1000.0 + np.testing.assert_allclose(trapezoid(basis, scaled, axis=0), 1.0, atol=1e-3) + + def test_ispline_monotonic_increasing(): X = np.linspace(0, 10, 200).reshape(-1, 1) transformer = ISplineTransformer(output_dim=8, include_bias=False)