R code for Bayes Shape-Restricted Regression Splines

BAYES SHAPE-RESTRICTED SPLINES


Non-parametric fits to a scatterplot with shape and smoothness assumptions.

If monotonicity is imposed, the quadratic B-spline basis functions are used. If the constraints include convexity (or concavity), the cubic B-spline basis functions are used.

The following code has arguments: type of shape constraint, x, y, k=#knots, and number of MCMC loops.

It returns the fit and point-wise confidence intervals at the x-values.

R code for Bayesian shape-restricted spline regression

Here is an example implementation, fitting a monotone convex curve to a dataset generated from an exponential function.

x=runif(50)

y=exp(2*x)+rnorm(50)

ans=brspl(5,x,y,3,5000)

plot(x,y)

lines(sort(x),ans$fit[order(x)])

lines(sort(x),ans$lower[order(x)],lty=3)

lines(sort(x),ans$upper[order(x)],lty=3)


Same as above but we allow for parametrically-modeled covariates to be included in the model.

For example, suppose we want to fit parallel increasing curves, so that a categorical covariate must be included in the smooth monotone regression. We construct k-1 dummy variables, where k is the number of levels of the covariate, and put these as columns in the matrix z. Note that the columns of z and the vector x must for a linearly independent set, and the one-vector must not be in the span of these vectors.

R code for Bayesian partial linear model with shape-restricted spline regression

Here is an example implementation, fitting three monotone convex curves to a data set generate from an exponential function.

x=1:50/50

x2=trunc(runif(50)*3)+1

z=matrix(0,nrow=50,ncol=2)

z[x2==2,1]=1

z[x2==3,2]=1

y=exp(2*x)+x2+rnorm(50)

ans=brspl(5,x,y,z,3,5000)

plot(x,y,col=x2)

lines(x[x2==1],ans$fit[x2==1])

lines(x[x2==2],ans$fit[x2==2],col=2)

lines(x[x2==3],ans$fit[x2==3],col=3)

meyer@stat.colostate.edu