R: Functions to Fit Univariate Break Point Regression Models
fitFixedKnots
R Documentation
Functions to Fit Univariate Break Point Regression Models
Description
There are two main functions here. The functions fitFixedKnots allows the fit a univariate regression using piecewise polynomials with "known" break points while the function fitFreeKnots estimates the break points.
the position of the interior knots for fitFixedKnots or starting values for fitFreeKnots
data
the data frame
degree
the degree if the piecewise polynomials
fixed
this is to be able to fit fixed break points
base
The basis for the piecewise polynomials, turn for truncated (default) and Bbase for B-base piecewise polynomials
trace
controlling the trace of of optim()
...
for extra arguments
Details
The functions fitFreeKnots() is loosely based on the curfit.free.knot() function of package
DierckxSpline of Sundar Dorai-Raj and Spencer Graves.
Value
The functions fitFixedKnots and fitFreeKnots return an object FixBreakPointsReg and
FreeBreakPointsReg respectively with the following items:
fitted.values
the fitted values of the model
residuals
the residuals of the model
df
the degrees of freedom fitted in the model
rss
the residuals sum of squares
knots
the knots used in creating the beta-function base
fixed
the fixed break points if any
breakPoints
the interior (estimated) break points (or knots)
coef
the coefficients of the linear part of the model
degree
the degree of the piecewise polynomial
y
the y variable
x
the x variable
w
the prior weights
Note
The prediction function in piecewise polynomials using the B-spline basis is tricky because by adding the newdata for x to the current one the B-basis function for the piecewise polynomials changes. This does not seems to be the case with the truncated basis, that is, when the option base="turn" is used (see the example below).
If the newdata are outside the range of the old x then there could a considerable discrepancies between the all fitted values and the predicted ones if the option base="Bbase" is used. The prediction function for the objects FixBreakPointsReg or FreeBreakPointsReg
has the option old.x.range=TRUE which allow the user two choices:
The first is to use the old end-points for the creation of the new B-basis which were determine from the original range of x. This choice is implemented as a default in the predict method for FixBreakPointsReg and FreeBreakPointsReg objects with the argument old.x.range=TRUE.
The second is to create new end-points from the new and old data x values. In this case
the range of x will be bigger that the original one if the newdata has values outside the original x range.
In this case (old.x.range=FALSE) the prediction could be possible better outside the x range but would not coincide with the original predictions i.e. fitted(model)
since basis have changed.
Dierckx, P. (1991) Curve and Surface Fitting with Splines, Oxford Science Publications
Stasinopoulos D. M. Rigby R.A. (2007) Generalized additive models for location scale and shape (GAMLSS) in R.
Journal of Statistical Software, Vol. 23, Issue 7, Dec 2007, http://www.jstatsoft.org/v23/i07.
Rigby, R. A. and Stasinopoulos D. M. (2005). Generalized additive models for location, scale and shape,(with discussion),
Appl. Statist., 54, part 3, pp 507-554.
Stasinopoulos D. M., Rigby R.A. and Akantziliotou C. (2006) Instructions on how to use the GAMLSS package in R.
Accompanying documentation in the current GAMLSS help files, (see also http://www.gamlss.org/)
Examples
# creating a linear + linear function
x <- seq(0,10, length.out=201)
knot <- 5
set.seed(12543)
mu <- ifelse(x<=knot,5+0.5*x,5+0.5*x+(x-knot))
y <- rNO(201, mu=mu, sigma=.5)
# plot the data
plot(y~x, xlim=c(-1,13), ylim=c(3,18))
# fit model using fixed break points
m1 <- fitFixedKnots(y, x, knots=5, degree=1)
knots(m1)
lines(fitted(m1)~x, col="red")
# now estimating the knot
m2 <- fitFreeKnots(y, x, knots=5, degree=1)
knots(m2)
summary(m2)
# now predicting
plot(y~x, xlim=c(-5,13), ylim=c(3,18))
lines(fitted(m2)~x, col="green", lwd=3)
points(-2:13,predict(m2, newdata=-2:13), col="red",pch = 21, bg="blue")
points(-2:13,predict(m2, newdata=-2:13, old.x.range=FALSE), col="red",pch = 21, bg="grey")
# fit different basis
m21 <- fitFreeKnots(y, x, knots=5, degree=1, base="Bbase")
deviance(m2)
deviance(m21) # should be identical
# predicting with m21
plot(y~x, xlim=c(-5,13), ylim=c(3,18))
lines(fitted(m21)~x, col="green", lwd=3)
points(-2:13,predict(m21, newdata=-2:13), col="red",pch = 21, bg="blue")
points(-2:13,predict(m21, newdata=-2:13, old.x.range=FALSE), col="red",pch = 21, bg="grey")