TermBuilder

Contents

TermBuilder#

class liesel_gam.TermBuilder(registry, prefix_names_by='', default_inference=MCMCSpec(<bound method IWLSKernel.untuned of <class 'liesel.goose.iwls.IWLSKernel'>>, self.kernel_group=None), default_scale_fn=(1.0, 0.005, 1.0))[source]#

Bases: object

Initializes structured additive model terms.

The terms returned by the methods of this class are all instances of liesel.model.Var, or of its subclasses.

Among other things, the term builder automatically assigns unique names to the created variables.

Parameters:
  • registry (PandasRegistry) – Provides an interface to a data frame used to set up the model terms.

  • prefix_names_by (str, default: '') – Names created by this TermBuilder will be prefixed by the string supplied here.

  • default_inference (Any | None, default: MCMCSpec(<bound method IWLSKernel.untuned of <class 'liesel.goose.iwls.IWLSKernel'>>, self.kernel_group=None)) – Defines the default inference specification for terms created by this builder. Note that this inference is only used for the coefficient variables of the terms created by this builder (StrctTerm.coef), not for the scale variables (StrctTerm.scale). The default gs.MCMCSpec(gs.IWLSKernel.untuned) is an iteratively-reweighted least squares kernel without step size tuning, see liesel.goose.IWLSKernel.untuned().

  • default_scale_fn (Callable[[], Var] | VarIGPrior, default: VarIGPrior(concentration=1.0, scale=0.005, value=1.0)) – A function or VarIGPrior object that defines the default scale for structured additive terms initialized by this builder. If this is a function, it must take no arguments and return a liesel.model.Var that acts as the scale. If it is a VarIGPrior, the default scale will be scale = sqrt(var), where var ~ InverseGamma(concentration, scale), with concentration and scale given by the VarIGPrior object. For most terms, this will mean that a fitting Gibbs sampler can be automatically set up for var. The exceptions to this rule are tf() and tx(). Note that, if you supply a custom default scale function, you should make sure that the inference attribute of your custom scale is defined, otherwise your custom scale may not be included in MCMC sampling. The default is VarIGPrior(1.0, 0.005), which leads to an inverse Gamma prior on the level of the variance parameter, i.e. \(\tau^2 \sim \operatorname{InverseGamma}(1.0, 0.005)\).

See also

BasisBuilder

Initializes Basis objects with penalty matrices.

Notes

The terms created by this builder generally have the form

\[s(\mathbf{x}_i) = \sum_{j=1}^J B_j(\mathbf{x}_i) \beta_j = \mathbf{b}(\mathbf{x}_i)^\top \boldsymbol{\beta}\]

where

  • \(i=1, \dots, N\) is the observation index,

  • \(\mathbf{x}_i^\top = [x_{i,1}, \dots, x_{i,M}]\) are covariate observations, where \(M\) denotes the number of covariates,

  • \(\mathbf{b}(\mathbf{x}_i)^\top = [B_1(\mathbf{x}_i), \dots, B_J(\mathbf{x}_i)]\) are a set of basis function evaluations, and

  • \(\boldsymbol{\beta}^\top = [\beta_1, \dots, \beta_J]\) are the corresponding coefficients.

In many cases, \(\mathbf{x}_i\) will consist of only one covariate, except for linear effects (lin(), slin()) or tensor product smooths (tf(), tx()).

The basis matrix for such a term is

\[\begin{split}\mathbf{B} = \begin{bmatrix} \mathbf{b}(\mathbf{x}_1)^\top \\ \vdots \\ \mathbf{b}(\mathbf{x}_N)^\top \end{bmatrix}.\end{split}\]

The coefficient receives a potentially rank-deficient multivariate normal prior

\[p(\boldsymbol{\beta}) \propto \left( \frac{1}{\tau^2}\right)^{\operatorname{rk}(\mathbf{K})/2} \exp \left( - \frac{1}{\tau^2} \boldsymbol{\beta}^\top \mathbf{K} \boldsymbol{\beta} \right)\]

with the potentially rank-deficient penalty matrix \(\mathbf{K}\) of rank \(\operatorname{rk}(\mathbf{K})\). The variance parameter \(\tau^2\) acts as an inverse smoothing parameter.

The choice of basis functions \(B_j\) and penalty matrix \(\mathbf{K}\) determines the nature of the term.

Sampling specification

By default, the coefficients for each term created with this termbuilder will be equipped with the TermBuilder’s default inference specification. This default specification corresponds to a blockwise sampling scheme, with each term’s coefficients forming a single block. See the examples for how to modify the default blockwise setup.

Overview of commonly used terms

Note

Basic terms

  • lin() : Linear term.

  • slin() : Linear term with iid penalty (ridge prior).

  • ps() : P-spline.

  • tp() : Thin plate spline.

  • ri() : Random intercept.

  • mrf() : Markov random field (discrete spatial effect).

  • kriging() : Low-rank gaussian process with fixed range.

Combined terms and tensor products

  • rs() : Random slope.

  • vc() : Varying coefficient.

  • tx() : Tensor product interaction without main effects.

  • tf() : Full tensor product with main effects.

Specialized smooths

  • np() : P-splines without linear trend.

  • cp() : Cyclic P-splines

Custom smooths

  • f() : Supply your own basis function and penalty matrix.

  • StrctTerm : Initialize a term independetly, potentially supplying a constant basis matrix and your own penalty matrix.

Tip

If your model somewhere contains a categorical variable, pay attention to the method labels_to_integers(); this helps you bring a newdata dictionary into a form understood by liesel.model.Model.predict() easily by turning string labels into their integer representations.

Examples

Basic example using defaults, initializing a P-spline:

>>> import liesel_gam as gam
>>> df = gam.demo_data(100)
>>> tb = gam.TermBuilder.from_df(df)
>>> tb.ps("x_nonlin", k=20)
StrctTerm(name="ps(x_nonlin)")

Changing the default inference. This still means, each term is sampled in a separate block.

>>> import liesel.goose as gs
>>> import liesel_gam as gam
>>> df = gam.demo_data(100)
>>> tb = gam.TermBuilder.from_df(df, default_inference=gs.MCMCSpec(gs.NUTSKernel))
>>> tb.ps("x_nonlin", k=20)
StrctTerm(name="ps(x_nonlin)")

Changing default inference, such that all terms’ coefficients are collected in the same block. Note that the scales are still sampled individually; by default with Gibbs kernels.

>>> import liesel.goose as gs
>>> import liesel_gam as gam
>>> df = gam.demo_data(100)
>>> tb = gam.TermBuilder.from_df(
...     df,
...     default_inference=gs.MCMCSpec(gs.NUTSKernel, kernel_group="1"),
... )
>>> tb.ps("x_nonlin", k=20)
StrctTerm(name="ps(x_nonlin)")

Changing the default scale. Here, we change the default scales to have a HalfNormal prior, and sample them on inverse-softplus level using independent IWLS kernels.

>>> import jax.numpy as jnp
>>> import liesel.model as lsl
>>> import tensorflow_probability.substrates.jax.bijectors as tfb
>>> import tensorflow_probability.substrates.jax.distributions as tfd
>>> import liesel_gam as gam
>>> def scale_fn():
...     prior = lsl.Dist(
...         tfd.HalfNormal,
...         scale=jnp.array(20.0),
...     )
...     scale = lsl.Var.new_param(
...         jnp.array(0.1),
...         distribution=prior,
...         name="{x}",  # placeholder for the automatically generated name
...     )
...     scale.transform(
...         tfb.Softplus(),
...         inference=gs.MCMCSpec(gs.IWLSKernel.untuned),  # inference for scales
...         name="h({x})",  # {x} is a placeholder
...     )
...     return scale
>>> df = gam.demo_data(100)
>>> tb = gam.TermBuilder.from_df(df, default_scale_fn=scale_fn)
>>> tb.ps("x_nonlin", k=20)
StrctTerm(name="ps(x_nonlin)")

Using a name prefix:

>>> import liesel_gam as gam
>>> df = gam.demo_data(100)
>>> tb = gam.TermBuilder.from_df(df, prefix_names_by="loc.")
>>> tb.ps("x_nonlin", k=20)
StrctTerm(name="loc.ps(loc.x_nonlin)")

If you don’t want the name prefix to appear on the covariate names, too, initialize the PandasRegistry individually. This way, you can for example use the same registry for two TermBuilder instances.

>>> import liesel_gam as gam
>>> df = gam.demo_data(100)
>>> registry = gam.PandasRegistry(df)
>>> tb1 = gam.TermBuilder(registry, prefix_names_by="loc.")
>>> tb1.ps("x_nonlin", k=20)
StrctTerm(name="loc.ps(x_nonlin)")
>>> tb2 = gam.TermBuilder(registry, prefix_names_by="scale.")
>>> tb2.ps("x_nonlin", k=20)
StrctTerm(name="scale.ps(x_nonlin)")

Methods

bs

B-spline with integrated squared derivative penalties.

cc

Cyclic version of cubic regression spline.

cp

Cyclic P-spline.

cr

Cubic regression spline.

cs

Cubic regression spline with additional null space shrinkage.

f

General StrctTerm, initialized by passing a custom basis function.

from_df

Initializes a TermBuilder from a pandas dataframe.

from_dict

Initializes a TermBuilder from a dictionary that holds the data.

init_scale

Initializes a scale variable with a term-related name.

kriging

Gaussian process models with a fixed range parameter in a basis-penalty-parameterization, often referred to as Kriging.

labels_to_integers

Processes a newdata dictionary, replacing labels of caterogical variables with their integer representation, such that they can be understood by liesel.model.Model.predict().

lin

Linear term.

mrf

Gaussian Markov random field.

np

P-spline without linear trend.

ps

P-spline: A B-spline with a discrete penalty matrix.

ri

Random intercept.

rs

Random slope.

slin

Linear term with an identity penalty matrix, leading to a ridge prior.

tf

General full anisotropic tensor product term, including main effects and interactions.

tp

Thin plate spline.

ts

Thin plate spline with additional null space penalty.

tx

General anisotropic tensor product interaction term without main effects.

vc

Varying coefficient term.