Spline features#

The modelling tools included in ISLP allow for construction of spline functions of features.

Force rebuild

import numpy as np
from ISLP import load_data
from ISLP.models import ModelSpec, ns, bs
Carseats = load_data('Carseats')
Carseats.columns
Index(['Sales', 'CompPrice', 'Income', 'Advertising', 'Population', 'Price',
       'ShelveLoc', 'Age', 'Education', 'Urban', 'US'],
      dtype='str')

Let’s make a term representing a cubic spline for Population. We’ll use knots based on the deciles.

knots = np.percentile(Carseats['Population'], np.linspace(10, 90, 9))
knots
array([ 58.9, 110.4, 160. , 218.6, 272. , 317.8, 366. , 412.2, 467. ])
bs_pop = bs('Population', internal_knots=knots, degree=3)

The object bs_pop does not refer to any data yet, it must be included in a ModelSpec object and fit using the fit method.

design = ModelSpec([bs_pop], intercept=False)
py_features = np.asarray(design.fit_transform(Carseats))

Compare to R#

We can compare our polynomials to a similar function in R

%load_ext rpy2.ipython
---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
Cell In[6], line 1
----> 1 get_ipython().run_line_magic('load_ext', 'rpy2.ipython')

File /usr/lib/python3.12/importlib/__init__.py:90, in import_module(name, package)
     88             break
     89         level += 1
---> 90 return _bootstrap._gcd_import(name[level:], package, level)

File <frozen importlib._bootstrap>:1387, in _gcd_import(name, package, level)

File <frozen importlib._bootstrap>:1360, in _find_and_load(name, import_)

File <frozen importlib._bootstrap>:1310, in _find_and_load_unlocked(name, import_)

File <frozen importlib._bootstrap>:488, in _call_with_frames_removed(f, *args, **kwds)

File <frozen importlib._bootstrap>:1387, in _gcd_import(name, package, level)

File <frozen importlib._bootstrap>:1360, in _find_and_load(name, import_)

File <frozen importlib._bootstrap>:1324, in _find_and_load_unlocked(name, import_)

ModuleNotFoundError: No module named 'rpy2'

We’ll recompute these features using bs in R. The default knot selection of the ISLP and R version are slightly different so we just fix the set of internal knots.

%%R -i Carseats,knots -o R_features
library(splines)
R_features = bs(Carseats$Population, knots=knots, degree=3)
UsageError: Cell magic `%%R` not found.
np.linalg.norm(py_features - R_features)
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[8], line 1
----> 1 np.linalg.norm(py_features - R_features)

NameError: name 'R_features' is not defined

Underlying model#

As for poly, the computation of the B-splines is done by a special sklearn transformer.

bs_pop
Feature(variables=('Population',), name='bs(Population, internal_knots=[ 58.9 110.4 160.  218.6 272.  317.8 366.  412.2 467. ], degree=3)', encoder=BSpline(internal_knots=array([ 58.9, 110.4, 160. , 218.6, 272. , 317.8, 366. , 412.2, 467. ]),
        lower_bound=np.float64(10.0), upper_bound=np.float64(509.0)), use_transform=True, pure_columns=False, override_encoder_colnames=True)

Natural splines#

Natural cubic splines are also implemented.

ns_pop = ns('Population', internal_knots=knots)
design = ModelSpec([ns_pop], intercept=False)
py_features = np.asarray(design.fit_transform(Carseats))
%%R -o R_features
library(splines)
R_features = ns(Carseats$Population, knots=knots)
UsageError: Cell magic `%%R` not found.
np.linalg.norm(py_features - R_features)
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
Cell In[12], line 1
----> 1 np.linalg.norm(py_features - R_features)

NameError: name 'R_features' is not defined

Intercept#

Looking at py_features we see it contains columns: [Population**i for i in range(1, 4)]. That is, it doesn’t contain an intercept, the order 0 term. This can be include with intercept=True. This means that the column space includes an intercept, though there is no specific column labeled as intercept.

bs_int = ns('Population', internal_knots=knots, intercept=True)
design = ModelSpec([bs_int], intercept=False)
py_int_features = np.asarray(design.fit_transform(Carseats))
py_int_features.shape, py_features.shape
((400, 11), (400, 10))