# FactorSmoothBasis


Factor-smooth interaction basis (per-level smooth with shared penalties).


Usage

``` python
FactorSmoothBasis(
    k=10,
    xt="tp",
    **marginal_kwargs,
)
```


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](RandomEffectBasis.md#whittaker.RandomEffectBasis) for level differences) when each level's mean *shape*, not merely its overall level, is expected to differ appreciably.


## Parameters


`k: int = ``10`  
Number 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: Any`  
Additional keyword arguments forwarded to the marginal basis constructor for the chosen `xt` (e.g. `m` for the spline order used by `"tp"` and `"ps"` bases).


## Notes

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:

1.  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.

2.  **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.


## Examples


``` python
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)


## Attributes

| Name | Description |
|----|----|
| [k](#k) | Marginal basis dimension per factor level. |
| [levels](#levels) | Sorted array of unique factor levels retained during `fit()`. |
| [marginal_basis](#marginal_basis) | The fitted marginal basis object shared, in form, across all levels. |
| [n_basis](#n_basis) | Total number of basis functions. |
| [n_levels](#n_levels) | Number of factor (grouping) levels seen during `fit()`. |

------------------------------------------------------------------------


### k


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.


------------------------------------------------------------------------


### levels


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.


------------------------------------------------------------------------


### marginal_basis


The fitted marginal basis object shared, in form, across all levels.


`marginal_basis: SmoothBasis`


This single fitted basis (e.g. a [TPRS](TPRS.md#whittaker.TPRS), [CRS](CRS.md#whittaker.CRS), or [PSpline](PSpline.md#whittaker.PSpline) instance) supplies the shape of the per-level smooths; it is evaluated once per level, each time multiplied by that level's indicator, to build the block-diagonal `basis_matrix()`.


------------------------------------------------------------------------


### n_basis


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.


#### Raises


`RuntimeError`  
If accessed before `fit()` has been called.


------------------------------------------------------------------------


### n_levels


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`.


## Methods

| Name | Description |
|----|----|
| [basis_matrix()](#basis_matrix) | Build the block-diagonal basis matrix. |
| [fit()](#fit) | Fit the factor-smooth basis. |
| [identifiability_constraints()](#identifiability_constraints) | Return `None`: no additional identifiability constraint is needed. |
| [null_space_dimension()](#null_space_dimension) | Return the dimension of the basis's unpenalized null space. |
| [penalty_matrices()](#penalty_matrices) | Return the full set of `1 + M` penalty matrices. |
| [penalty_matrix()](#penalty_matrix) | Return only the shared wiggliness penalty (first entry of `penalty_matrices()`). |

------------------------------------------------------------------------


### basis_matrix()


Build the block-diagonal basis matrix.


Usage

``` python
basis_matrix(
    x_numeric,
    factor=None,
    **kwargs,
)
```


#### Parameters


`x_numeric: NDArray`  
Numeric covariate, shape `(n,)`.

`factor: NDArray | None = None`  
Factor variable, shape `(n,)`.


#### Returns


`NDArray`  
Design matrix of shape `(n, n_levels * k_marginal)`.


------------------------------------------------------------------------


### fit()


Fit the factor-smooth basis.


Usage

``` python
fit(
    x_numeric,
    factor=None,
    **kwargs,
)
```


#### Parameters


`x_numeric: NDArray`  
Numeric covariate, shape `(n,)`.

`factor: NDArray | None = None`  
Factor (grouping) variable, shape `(n,)`. Can be strings or integers.


------------------------------------------------------------------------


### identifiability_constraints()


Return `None`: no additional identifiability constraint is needed.


Usage

``` python
identifiability_constraints()
```


Unlike [RandomEffectBasis](RandomEffectBasis.md#whittaker.RandomEffectBasis) or [MRFBasis](MRFBasis.md#whittaker.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.


------------------------------------------------------------------------


### null_space_dimension()


Return the dimension of the basis's unpenalized null space.


Usage

``` python
null_space_dimension()
```


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`.


#### Returns


`int`  
Always `0`.


------------------------------------------------------------------------


### penalty_matrices()


Return the full set of `1 + M` penalty matrices.


Usage

``` python
penalty_matrices()
```


#### Returns


`list[NDArray]`  
A list of `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.


------------------------------------------------------------------------


### penalty_matrix()


Return only the shared wiggliness penalty (first entry of `penalty_matrices()`).


Usage

``` python
penalty_matrix()
```


This exists for compatibility with the single-penalty [SmoothBasis](SmoothBasis.md#whittaker.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.


#### Returns


`NDArray`  
The `(k_total, k_total)` wiggliness penalty `I_L ⊗ S_marginal`.