gamlss2

Monotonic P-Splines for GAMLSS

Description

Constructs a monotonic penalized B-spline (P-spline) model term for use with gamlss2. The smooth can be constrained to be monotonically increasing or decreasing. The smoothing parameter is either supplied by the user or selected automatically during estimation.

Usage

ms(x, mono = c("up", "down"), k = 10,
  criterion = "aicc", lambda = NULL,
  lambda.min = 1e-8, lambda.max = 1e8,
  lambda.tol = 1e-4, ...)

Arguments

x

A numeric covariate to be represented by a monotonic P-spline.

mono

Character specifying the direction of monotonicity. Use “up” for a monotonically increasing effect or “down” for a monotonically decreasing effect.

k

Integer giving the dimension of the B-spline basis. The default is 10.

criterion

Character specifying the criterion used for automatic smoothing parameter selection. Available criteria are “aicc”, “aic”, “gaic”, “bic”, and “gcv”. The default is “aicc”.

lambda

Optional positive smoothing parameter. If NULL, the smoothing parameter is selected automatically. If supplied, lambda is kept fixed during fitting.

lambda.min, lambda.max

Positive lower and upper limits for automatic smoothing parameter selection. The search is performed on the logarithmic scale.

lambda.tol

Tolerance used in the one-dimensional optimization of the smoothing parameter.

…

Further arguments passed to s when constructing the P-spline basis.

Details

Function ms() defines a special model term for gamlss2. The underlying basis is an mgcv P-spline basis constructed with bs = “ps”.

For an increasing smooth, adjacent B-spline coefficients are constrained to satisfy

\(\beta_{r+1} - \beta_r \ge 0,\)

whereas for a decreasing smooth the inequalities are reversed. Because B-spline basis functions are non-negative and form a partition of unity, ordering the spline coefficients in this way produces a monotonic fitted effect.

Estimation is based on penalized weighted least squares. For a given smoothing parameter \(\lambda\), the coefficient vector minimizes

(z - X)^W (z - X) + ^S

subject to the monotonicity constraints. The constrained problem is solved with an active-set algorithm. Weighted cross-products are computed once per backfitting update, so smoothing parameter selection only requires operations on matrices having dimension determined by k.

If lambda = NULL, the smoothing parameter is selected by a one-dimensional search on log(lambda) using the criterion specified by criterion.

The fitted term is centered before it is returned to the additive backfitting algorithm.

Value

The constructor ms() returns an object of classes “special” and “ms” containing the P-spline design and penalty matrices and information required for fitting and prediction.

After fitting with gamlss2, the corresponding fitted special term contains the following main components:

coefficients: estimated B-spline coefficients.
fitted.values: centered fitted contribution of the term.
edf: effective degrees of freedom.
lambdas: selected or fixed smoothing parameter.
vcov: covariance/influence matrix conditional on the active constraints.
df: residual degrees of freedom for the term, n - edf.

Additional components are stored for prediction and for warm starts in subsequent backfitting iterations.

References

Eilers, P. H. C. and Marx, B. D. (1996). “Flexible Smoothing with B-Splines and Penalties.” Statistical Science, 11(2), 89–121. doi:10.1214/ss/1038425655

Bollaerts, K., Eilers, P. H. C. and van Mechelen, I. (2006). “Simple and Multiple P-splines Regression with Shape Constraints.” British Journal of Mathematical and Statistical Psychology, 59(2), 451–469. doi:10.1348/00071105X84293

Examples

library("gamlss2")

air <- subset(airquality, !is.na(Ozone))
air <- air[order(air$Temp), ]

## monotone Temp effect for mu and sigma
## Wind is modeled by an unconstrained smooth
f <- Ozone ~ s(Wind) + ms(Temp, "up") | ms(Temp, "down")

## estimate model
m <- gamlss2(f, data = air, family = GA)

## model summary
summary(m)

## visualize estimated effects
plot(m)