Generates the nonnegative natural cubic spline basis matrix, the
corresponding integrals (from the left boundary knot), or derivatives of
given order. Each basis is assumed to follow a linear trend for x
outside of boundary.
Usage
naturalSpline(
x,
df = NULL,
knots = NULL,
intercept = FALSE,
Boundary.knots = NULL,
derivs = 0L,
integral = FALSE,
...
)
Arguments
- x
The predictor variable. Missing values are allowed and will be returned as they are.
- df
Degree of freedom that equals to the column number of returned matrix. One can specify
df
rather thanknots
, then the function choosesdf - 1 - as.integer(intercept)
internal knots at suitable quantiles ofx
ignoring missing values and thosex
outside of the boundary. Thus,df
must be greater than or equal to2
. If internal knots are specified viaknots
, the specifieddf
will be ignored.- knots
The internal breakpoints that define the splines. The default is
NULL
, which results in a basis for ordinary polynomial regression. Typical values are the mean or median for one knot, quantiles for more knots.- intercept
If
TRUE
, the complete basis matrix will be returned. Otherwise, the first basis will be excluded from the output.- Boundary.knots
Boundary points at which to anchor the splines. By default, they are the range of
x
excludingNA
. If bothknots
andBoundary.knots
are supplied, the basis parameters do not depend onx
. Data can extend beyondBoundary.knots
.- derivs
A nonnegative integer specifying the order of derivatives of natural splines. The default value is
0L
for the spline basis functions.- integral
A logical value. The default value is
FALSE
. IfTRUE
, this function will return the integrated natural splines from the left boundary knot.- ...
Optional arguments that are not used.
Value
A numeric matrix of length(x)
rows and df
columns if df
is specified or length(knots) + 1 +
as.integer(intercept)
columns if knots
are specified instead.
Attributes that correspond to the arguments specified are returned for
usage of other functions in this package.
Details
It is an implementation of the natural spline basis based on B-spline basis, which utilizes the close-form null space that can be derived from the recursive formula for the second derivatives of B-splines. The constructed spline basis functions are intended to be nonnegative within boundary with second derivatives being zeros at boundary knots.
A similar implementation is provided by splines::ns
, which uses QR
decomposition to find the null space of the second derivatives of B-spline
basis at boundary knots. However, there is no guarantee that the resulting
basis functions are nonnegative within boundary.
Examples
library(splines2)
x <- seq.int(0, 1, 0.01)
knots <- c(0.3, 0.5, 0.6)
## natural spline basis
nsMat0 <- naturalSpline(x, knots = knots, intercept = TRUE)
## integrals
nsMat1 <- naturalSpline(x, knots = knots, intercept = TRUE, integral = TRUE)
## first derivatives
nsMat2 <- naturalSpline(x, knots = knots, intercept = TRUE, derivs = 1)
## second derivatives
nsMat3 <- naturalSpline(x, knots = knots, intercept = TRUE, derivs = 2)
op <- par(mfrow = c(2, 2), mar = c(2.5, 2.5, 0.2, 0.1), mgp = c(1.5, 0.5, 0))
matplot(x, nsMat0, type = "l", ylab = "basis")
matplot(x, nsMat1, type = "l", ylab = "integral")
matplot(x, nsMat2, type = "l", ylab = "1st derivative")
matplot(x, nsMat3, type = "l", ylab = "2nd derivative")
par(op) # reset to previous plotting settings
## use the deriv method
all.equal(nsMat0, deriv(nsMat1), check.attributes = FALSE)
#> [1] TRUE
all.equal(nsMat2, deriv(nsMat0))
#> [1] TRUE
all.equal(nsMat3, deriv(nsMat2))
#> [1] TRUE
all.equal(nsMat3, deriv(nsMat0, 2))
#> [1] TRUE