Factor-smooth interaction basis (per-level smooth with shared penalties).
Usage
A factor-smooth interaction fits a separate curve of a numeric covariate for every level of a grouping factor — one trajectory per subject in a longitudinal study, one seasonal curve per site, one dose-response curve per batch — while pooling information across levels by sharing the same wiggliness penalty (and, optionally, coefficient-level shrinkage) across all of them. It is equivalent to mgcv’s bs="fs" basis and is the natural GAM analogue of a random-slope mixed model. Choose this basis over fitting k independent smooths (one per level) when the levels are numerous, some levels have little data, and it is desirable to borrow strength across levels via a shared smoothing parameter; choose it over a single shared smooth (with an additional RandomEffectBasis for level differences) when each level’s mean shape, not merely its overall level, is expected to differ appreciably.
k: int = 10Number of basis functions for the marginal smooth per level. Applies to every level identically — all per-level smooths share the same basis dimension. Larger k allows more wiggly per-level curves at the cost of more coefficients (n_levels * k total); since the smoothing parameter (not k) ultimately controls how wiggly the fitted curves are, k mainly needs to be large enough not to unduly constrain the shape. The default is 10.
xt: str = "tp"Marginal basis type used for each level’s smooth. One of "tp" (thin plate regression spline, a good general-purpose default), "cr" (cubic regression spline, cheaper for a single covariate with many knots), or "ps" (P-spline, useful when a difference penalty on B-spline coefficients is preferred). The default is "tp".
**marginal_kwargs: Anyxt (e.g. m for the spline order used by "tp" and "ps" bases).
Let there be L factor levels and let the fitted marginal basis (shared in form, but evaluated per level) have k_m basis functions with null-space dimension M. The full basis matrix is block-diagonal by level:
\mathbf{B} = \begin{bmatrix} \mathbf{1}_{[\text{level}=1]} \odot \mathbf{B}_m & \mathbf{1}_{[\text{level}=2]} \odot \mathbf{B}_m & \cdots & \mathbf{1}_{[\text{level}=L]} \odot \mathbf{B}_m \end{bmatrix},
where B_m is the marginal basis matrix evaluated at the numeric covariate and \mathbf{1}_{[\text{level}=\ell]} is the indicator for observations belonging to level \ell (each row of B is nonzero only in the block for its own level); the total number of columns is L * k_m. The penalty structure has 1 + M components:
A shared wiggliness penalty, replicated identically across all levels via the Kronecker structure
\mathbf{S}_{\text{wiggle}} = \mathbf{I}_L \otimes \mathbf{S}_m ,
where S_m is the marginal basis’s own penalty matrix. A single smoothing parameter controls how wiggly every level’s curve is allowed to be.
One penalty per marginal null-space component (there are M of them, e.g. M=2 for a thin plate spline with m=2 — the constant and linear components). For each null-space eigenvector v of S_m, the corresponding penalty places \operatorname{outer}(v, v) in every level’s diagonal block, so that this penalty shrinks that low-order component of the curve (e.g. each level’s intercept, or each level’s linear trend) toward a common value across levels — a random-intercept/random-slope penalty for exactly the marginal basis’s unpenalized directions.
Because every basis coefficient is touched by at least one of these 1 + M penalties (the wiggliness penalty covers the range space and the null-space penalties cover what the wiggliness penalty leaves unpenalized), null_space_dimension() is always 0 and no identifiability_constraints() are required — the basis is fully penalized and does not collide with a fixed intercept. fit() raises ValueError if fewer than 2 factor levels are present, since factor-smooth interactions require differentiating between levels.
import numpy as np
from whittaker.smooths.factor_smooth import FactorSmoothBasis
rng = np.random.default_rng(0)
n = 100
subject = rng.choice(["s1", "s2", "s3"], size=n)
x_numeric = rng.uniform(0, 1, n)
basis = FactorSmoothBasis(k=8).fit(x_numeric, subject)
B = basis.basis_matrix(x_numeric, subject)
penalties = basis.penalty_matrices()
B.shape, len(penalties)((100, 24), 3)| Name | Description |
|---|---|
| k | Marginal basis dimension per factor level. |
| levels |
Sorted array of unique factor levels retained during fit().
|
| marginal_basis | The fitted marginal basis object shared, in form, across all levels. |
| n_basis | Total number of basis functions. |
| n_levels |
Number of factor (grouping) levels seen during fit().
|
Marginal basis dimension per factor level.
k: int
Before fit() is called, returns the requested value passed to __init__ (which may be -1 to defer to the marginal basis’s own default). After fit(), returns the marginal basis’s actual, fitted dimension.
Sorted array of unique factor levels retained during fit().
levels: NDArray
Each level corresponds to one diagonal block of the basis matrix and one column range of k_marginal coefficients.
The fitted marginal basis object shared, in form, across all levels.
marginal_basis: SmoothBasis
Total number of basis functions.
n_basis: int
Equal to n_levels * k_marginal: the marginal basis dimension repeated once per factor level, since the full basis matrix is block-diagonal by level.
RuntimeErrorfit() has been called.
Number of factor (grouping) levels seen during fit().
n_levels: int
Determines the number of diagonal blocks in the basis matrix and, together with k_marginal, the total basis dimension n_basis.
| Name | Description |
|---|---|
| basis_matrix() | Build the block-diagonal basis matrix. |
| fit() | Fit the factor-smooth basis. |
| identifiability_constraints() |
Return None: no additional identifiability constraint is needed.
|
| null_space_dimension() | Return the dimension of the basis’s unpenalized null space. |
| penalty_matrices() |
Return the full set of 1 + M penalty matrices.
|
| penalty_matrix() |
Return only the shared wiggliness penalty (first entry of penalty_matrices()).
|
Build the block-diagonal basis matrix.
Usage
x_numeric: NDArrayNumeric covariate, shape (n,).
factor: NDArray | None = None(n,).
NDArray(n, n_levels * k_marginal).
Fit the factor-smooth basis.
Usage
x_numeric: NDArrayNumeric covariate, shape (n,).
factor: NDArray | None = None(n,). Can be strings or integers.
Return None: no additional identifiability constraint is needed.
Usage
Unlike RandomEffectBasis or MRFBasis, the null-space penalties here already shrink each level’s low-order components toward a common value, so there is no separate unpenalized constant that collides with a model intercept.
Return the dimension of the basis’s unpenalized null space.
Usage
The factor-smooth basis has no unpenalized null space: the shared wiggliness penalty I_L ⊗ S_marginal covers the range space of the marginal penalty in every level, and the M additional null-space (random-intercept/slope) penalties returned by penalty_matrices() cover exactly the directions the wiggliness penalty leaves unpenalized. Since every basis coefficient is touched by at least one penalty, this is always 0.
int0.
Return the full set of 1 + M penalty matrices.
Usage
list[NDArray]k_total x k_total matrices (k_total = n_levels * k_marginal): the first is the shared wiggliness penalty I_L ⊗ S_marginal, and the remaining M (the marginal basis’s null-space dimension) are random-effect-style penalties, one per null-space component of the marginal penalty, each shrinking that component’s coefficients toward a common value across levels.
Return only the shared wiggliness penalty (first entry of penalty_matrices()).
Usage
This exists for compatibility with the single-penalty SmoothBasis interface. It omits the M null-space (random-intercept/slope) penalties; use penalty_matrices() to access the full penalty structure needed for a proper fit.
NDArray(k_total, k_total) wiggliness penalty I_L ⊗ S_marginal.