mgcv (version 0.6-1)

mono.con: Monotonicity constraints for a cubic regression spline.

Description

Finds linear constraints sufficient for monotonicity (and optionally upper and/or lower boundedness) of a cubic regression spline. The basis representation assumed is that given by the gam, "cr" basis: that is the spline has a set of knots, which have fixed x values, but the y values of which constitute the parameters of the spline.

Usage

mono.con(x,up=TRUE,lower=NA,upper=NA)

Arguments

x

The array of knot locations.

up

If TRUE then the constraints imply increase, if FALSE then decrease.

lower

This specifies the lower bound on the spline unless it is NA in which case no lower bound is imposed.

upper

This specifies the upper bound on the spline unless it is NA in which case no upper bound is imposed.

Value

• The function returns a list containing constraint matrix A and constraint vector b.

Details

Consider the natural cubic spline passing through the points: ${x_i,p_i:i=1 \ldots n }$. Then it is possible to find a realtively small set of linear constraints on $\bf p$ sufficient to ensure monotonicity (and bounds if required): ${\bf Ap}\ge{\bf b}$. Details are given in Wood (1994). This function returns a list containing A and b.

References

Gill, P.E., Murray, W. and Wright, M.H. (1981) Practical Optimization. Academic Press, London.

Wood, S.N. (1994) Monotonic smoothing splines fitted by cross validation SIAM Journal on Scientific Computing 15(5):1126-1133

http://www.ruwpa.st-and.ac.uk/simon.html

See Also

mgcv pcls

Examples

# Fit a monotonic penalized regression spline .....

# Generate data from a monotonic truth.
set.seed(10);x<-runif(100)*4-1;x<-sort(x);
f<-exp(4*x)/(1+exp(4*x));y<-f+rnorm(100)*0.1;plot(x,y)
# Show regular spline fit (and save fitted object)
f.ug<-gam(y~s(x,k=10,bs="cr"));lines(x,fitted(f.ug))
# Create Design matrix, constriants etc. for monotonic spline....
gam.setup(y~s(x,k=10,bs="cr")-1)->G;GAMsetup(G)->G;F<-mono.con(G$xp);
G$Ain<-F$A;G$bin<-F$b;G$C<-matrix(0,0,0);G$sp<-f.ug$sp;
G$p<-G$xp;G$y<-y;G$w<-y*0+1

p<-pcls(G);  # fit spline (using s.p. from unconstrained fit)

# now modify the gam object from unconstrained fit a little, to use it
# for predicting and plotting constrained fit. 
p<-c(0,p);f.ug$coefficients<-p; 
lines(x,predict.gam(f.ug,newdata=data.frame(x=x)),col=2)