grpnet (version 0.5)

rk: Reproducing Kernel Basis

Description

Generate a reproducing kernel basis matrix for a nominal, ordinal, or polynomial smoothing spline.

Usage

rk(x, df = NULL, knots = NULL, m = NULL, intercept = FALSE, 
   Boundary.knots = NULL, warn.outside = TRUE, 
   periodic = FALSE, xlev = levels(x))

Value

Returns a basis function matrix of dimension n by df (plus 1 if an intercept is included) with the following attributes:

df

degrees of freedom

knots

knots for spline basis

m

derivative penalty order

intercept

was an intercept included?

Boundary.knots

boundary points of x

periodic

is the basis periodic?

xlev

factor levels (if applicable)

Arguments

x

the predictor vector of length n. Can be a factor, integer, or numeric, see Note.

df

the degrees of freedom, i.e., number of knots to place at quantiles of x. Defaults to 5 but ignored if knots are provided.

knots

the breakpoints (knots) defining the spline. If knots are provided, the df is defined as length(unique(c(knots, Boundary.knots))).

m

the derivative penalty order: 0 = ordinal spline, 1 = linear spline, 2 = cubic spline, 3 = quintic spline

intercept

should an intercept be included in the basis?

Boundary.knots

the boundary points for spline basis. Defaults to range(x).

warn.outside

if TRUE, a warning is provided when x values are outside of the Boundary.knots

periodic

should the spline basis functions be constrained to be periodic with respect to the Boundary.knots?

xlev

levels of x (only applicable if x is a factor)

Author

Nathaniel E. Helwig <[e-mail]>

Details

Given a vector of function realizations \(f\), suppose that \(f = X \beta\), where \(X\) is the (unregularized) spline basis and \(\beta\) is the coefficient vector. Let \(Q\) denote the postive semi-definite penalty matrix, such that \(\beta^\top Q \beta\) defines the roughness penalty for the spline. See Helwig (2017) for the form of \(X\) and \(Q\) for the various types of splines.

Consider the spectral parameterization of the form \(f = Z \alpha\) where $$Z = X Q^{-1/2}$$ is the regularized spline basis (that is returned by this function), and \(\alpha = Q^{1/2} \beta\) are the reparameterized coefficients. Note that \(X \beta = Z \alpha\) and \(\beta^\top Q \beta = \alpha^\top \alpha\), so the spectral parameterization absorbs the penalty into the coefficients (see Helwig, 2021, 2024).

Syntax of this function is designed to mimic the syntax of the bs function.

References

Helwig, N. E. (2017). Regression with ordered predictors via ordinal smoothing splines. Frontiers in Applied Mathematics and Statistics, 3(15), 1-13. tools:::Rd_expr_doi("10.3389/fams.2017.00015")

Helwig, N. E. (2021). Spectrally sparse nonparametric regression via elastic net regularized smoothers. Journal of Computational and Graphical Statistics, 30(1), 182-191. tools:::Rd_expr_doi("10.1080/10618600.2020.1806855")

Helwig, N. E. (2024). Precise tensor product smoothing via spectral splines. Stats, 7(1), 34-53. tools:::Rd_expr_doi("10.3390/stats7010003")

Examples

######***######   NOMINAL SPLINE BASIS   ######***######

x <- as.factor(LETTERS[1:5])
basis <- rk(x)
plot(1:5, basis[,1], t = "l", ylim = extendrange(basis))
for(j in 2:ncol(basis)){
  lines(1:5, basis[,j], col = j)
}


######***######   ORDINAL SPLINE BASIS   ######***######

x <- as.ordered(LETTERS[1:5])
basis <- rk(x)
plot(1:5, basis[,1], t = "l", ylim = extendrange(basis))
for(j in 2:ncol(basis)){
  lines(1:5, basis[,j], col = j)
}


######***######   LINEAR SPLINE BASIS   ######***######

x <- seq(0, 1, length.out = 101)
basis <- rk(x, m = 1)
plot(x, basis[,1], t = "l", ylim = extendrange(basis))
for(j in 2:ncol(basis)){
  lines(x, basis[,j], col = j)
}


######***######   CUBIC SPLINE BASIS   ######***######

x <- seq(0, 1, length.out = 101)
basis <- rk(x)
basis <- scale(basis)  # for visualization only!
plot(x, basis[,1], t = "l", ylim = extendrange(basis))
for(j in 2:ncol(basis)){
  lines(x, basis[,j], col = j)
}


######***######   QUINTIC SPLINE BASIS   ######***######

x <- seq(0, 1, length.out = 101)
basis <- rk(x, m = 3)
basis <- scale(basis)  # for visualization only!
plot(x, basis[,1], t = "l", ylim = extendrange(basis))
for(j in 2:ncol(basis)){
  lines(x, basis[,j], col = j)
}