Require that the spline be convex

Question

I need to fit the spline to the dataset, and the resulting function should be monotonically decreasing and convex. The data that I pass to splinefun is guaranteed to have these properties, but this does not guarantee that the resulting function is convex. Is there a way to fit a spline to a dataset and require the resulting function to be convex?

Answer 1

First, provide some example data:

x = c(0,1,2,3,4,5,6) y = c(2,1, 0.59, 0.27, 0.25, -0.23, -0.45) dat <- data.frame(x=x,y=y) 

Monotone Spline

We can make a monotone spline with splinefun(x,y,"monoH.FC") , as @fang suggests.

 # Setting up Monotonic Spline MonoSpline = splinefun(x,y,"monoH.FC") #Getting Ready for plotting Monotonic Spline xArray = seq(0,6,0.01) MonoResult = MonoSpline(xArray) 

Monotonic convex spline (with scam package)

For a monotonous convex spline, you need to use the scam package. Then we can:

 # Setting up Monotonic Convex Spline # install.packages("scam") require(scam) MonoConvexSpline <- scam(y~s(x,k=4,bs="mdcx",m=1),data=dat) MonoConvexSplinePredict =function(Test){ predict.scam(MonoConvexSpline,data.frame(x = Test)) } #Getting Ready for plotting Monotonic Convex Spline MonoConvexSplineResult = MonoConvexSplinePredict(xArray) 

Please note the following:

• The bs="mdcx" means that we want to reduce the convex spline. If you want to increase convex, decreasing concave, etc., then find the corresponding bs here .
• If you put non-montonal data in the splinefun(x,y,"monoH.FC") , we get an error.
• If you put non-convex data in the scam function, you will still get a spline. Data changes so that small convexes change to convex. There are no warnings about this, so be careful, because your data may look completely different. As an example, the graph below was executed using the same code as above, except for the decreasing concave function bs="mdcv" :

Monotonous convex spline (with cobs package)

 # Convex Cobs Spline library(cobs) spCobs = cobs(x , y, constraint = c("decrease", "convex"), nknots = 8) spCobsResults = predict(spCobs, xArray)[,2] 

Building everything

Then build them using

 Plot = qplot(xlab = "x", ylab = "y") Plot = Plot + geom_line(aes(xArray,MonoResult , colour = "Monotonic Spline" )) Plot = Plot + geom_line(aes(xArray,MonoConvexSplineResult, colour = "Monotonic Convex scam Spline")) Plot = Plot + geom_line(aes(xArray,spCobsResults , colour = "Monotonic Convex cobs Spline")) Plot 

Speed

• Predicting points using the scam function takes a lot more time than using a monotonic spline or cob spline. You can see it with the help of micro-business below.
• Alloys of splines and spline fraudsters take much more time than initially calculated than a common monotonous spline.

.

 # Prediction library(microbenchmark) microbenchmark( MonoSpline(xArray), predict.scam(MonoConvexSpline,data.frame(x = xArray)), predict(spCobs, xArray)[,2] ) Unit: microseconds expr min lq mean median uq max neval MonoSpline(xArray) 141.540 147.8175 223.3695 156.9490 167.9830 1593.456 100 predict.scam(MonoConvexSpline, data.frame(x = xArray)) 2778.655 2838.0095 3161.2282 2914.8665 3153.4285 6168.741 100 predict(spCobs, xArray)[, 2] 125.179 133.1690 155.1226 145.1535 162.2755 366.784 100 # Calculating Spline library(microbenchmark) microbenchmark( splinefun(x,y,"monoH.FC"), scam(y~s(x,k=4,bs="mdcx",m=1),data=dat), cobs(x , y, constraint = c("decrease", "convex"), nknots = 8) ) Unit: microseconds expr min lq mean median uq max neval splinefun(x, y, "monoH.FC") 90.175 127.462 411.6407 153.7155 198.993 24877.47 100 scam(y ~ s(x, k = 4, bs = "mdcx", m = 1), data = dat) 166769.270 196719.139 231631.5321 224372.7940 265074.525 355734.37 100 cobs(x, y, constraint = c("decrease", "convex"), nknots = 8) 145511.335 172887.618 203786.0940 202997.4795 228688.607 347661.29 100 

Answer 2

My other answer to this question showed a monotonous spline shape, as well as curls and curls. The problem with these shape-limited splines is that they are rather slow and do not necessarily interpolate all data points.

Spline Schumaker

I released a package that implements the Sumerak spline, which is monotone and convex / concave if the data is monotone and convex / concave. It is fast and interpolates all data points.

Example:

 #install.packages("schumaker") library(schumaker) x = seq(1,10) y = -log(x) xarray = seq(1,10,0.01) SchumSpline = schumaker::Schumaker(x,y) Schum0 = SchumSpline$Spline(xarray) Schum1 = SchumSpline$DerivativeSpline(xarray) Schum2 = SchumSpline$SecondDerivativeSpline(xarray) plot(xarray, Schum0, type = "l", col = 4, ylim = c(-3,1), main = "Schumaker Spline and first two derivatives", ylab = "Spline and derivatives", xlab = "x") points(x,y) lines(xarray, Schum1, col = 2) lines(xarray, Schum2, col = 3) abline(h = 0, col = 1) text(x=rep(8,8,8), y=c(-2, -0.5,+0.2), pos=4, labels=c('Spline', 'First Derivative', 'Second Derivative')) 

You can see that the second derivative is always positive (this is not true for a regular monotone spline. See package vignette).

Note that the spline will only be globally convex / concave if the data is globally convex / concave. This is unavoidable since it is an interpolating spline.

This spline is faster than cobs and fraudulent splines. This is slower to create than a monotonous spline, but faster to evaluate. Full speed tests can be found in the vignette.