Shape Constraints with BsplineQuantReg

Introduction

This vignette covers shape constraints in B-spline quantile regression using the Karlin-Studden SOCP formulation. Shape constraints allow you to incorporate prior knowledge about the function’s behavior (monotonicity, convexity, etc.) into the estimation process.

Shape constraint for quantile regression with cubic splines is initially the main purpose of this package, based on the paper of Abbes (2025). It has been extended to other degrees, and the package also allows to switch to ordinary mean square regression.

Available Demos

The package includes several demos that illustrate shape constraints:

# Example of available demos
#dev.new()
# Run specific demos:
# demo(package = "BsplineQuantReg","monotonicity")      #- Increasing/decreasing constraints
# demo(package = "BsplineQuantReg","convexity")         #- Convexity/concavity constraints  
# demo(package = "BsplineQuantReg","der3")              #- Third derivative constraints
# demo(package = "BsplineQuantReg","comprehensive")     #- All constraint types
# demo(package = "BsplineQuantReg","temperature") #and
#demo(package = "BsplineQuantReg",temperature2) #- Real-world example with temperature data

Uniform Constraints

Monotonicity Constraints

Monotonicity constraints force the fitted function to be non-decreasing (monot = 1) or non-increasing (monot = -1).

# Generate data
set.seed(42)
n <- 100
x <- seq(0, 1, length.out = n)
y <- 2*x + 0.5*sin(4*pi*x) + 0.1*rnorm(n)
knots <- quantile(x, probs = seq(0, 1, length.out = 10))

# Unconstrained
fit_uncon <- quantile_spline(x, y, knots, tau = 0.5)

# Increasing constraint
fit_inc <- quantile_spline(x, y, knots, tau = 0.5, monot = 1)

# Decreasing constraint
fit_dec <- quantile_spline(x, y, knots, tau = 0.5, monot = -1)

# Compare
x_eval <- seq(0, 1, length.out = 200)
plot(x, y, pch = 16, cex = 0.5, col = "gray", main = "Monotonicity Constraints")
lines(x_eval, spline_eval(fit_uncon, x_eval), col = "red", lwd = 2, lty = 2)
lines(x_eval, spline_eval(fit_inc, x_eval), col = "blue", lwd = 2)
lines(x_eval, spline_eval(fit_dec, x_eval), col = "green", lwd = 2)
legend("topleft", legend = c("Unconstrained", "Increasing", "Decreasing"),
       col = c("red", "blue", "green"), lty = c(2, 1, 1), lwd = 2)

Convexity Constraints

Convexity constraints force the second derivative to be non-negative (convcons = 1) or non-positive (convcons = -1).

# Convex function
y_conv <- (x - 0.5)^2 + 0.05*rnorm(n)

fit_convex <- quantile_spline(x, y_conv, knots, tau = 0.5, convcons = 1)

plot(x, y_conv, pch = 16, cex = 0.5, col = "gray", main = "Convexity Constraint")
lines(x_eval, spline_eval(fit_convex, x_eval), col = "blue", lwd = 2)

Third Derivative Constraints

Third derivative constraints are not commonly intuitive. They control the variations of the second derivative, that is to say the evolution of curvature. Some examples show that the control of its sign allows sometimes enforcing smoothness of the curvature, and strongly enforces smoothness of the curve.

# Cubic function with varying third derivative
y_cubic <- x^3 - 1.5*x^2 + 0.5*x + 0.05*rnorm(n)

fit_der3_pos <- quantile_spline(x, y_cubic, knots, tau = 0.5, der3cons = 1)

plot(x, y_cubic, pch = 16, cex = 0.5, col = "gray", main = "Third Derivative Constraint")
lines(x_eval, spline_eval(fit_der3_pos, x_eval), col = "blue", lwd = 2)

Partial Constraints

Partial constraints apply monotonicity or convexity only on specific intervals, not the entire domain.

Partial Monotonicity

# Increasing only on first half
monot_partial <- c(rep(1, 5), rep(0, 4))  # 5 intervals increasing, 4 unconstrained

fit_partial <- quantile_spline(x, y, knots, tau = 0.5, monot = monot_partial)

plot(x, y, pch = 16, cex = 0.5, col = "gray", main = "Partial Monotonicity")
lines(x_eval, spline_eval(fit_partial, x_eval), col = "blue", lwd = 2)
abline(v = knots[6], col = "red", lty = 2)  # Boundary of constrained region

Partial Convexity

# Convex on right half only
conv_partial <- rep(0, length(knots))
for (i in 1:length(knots)) {
  if (knots[i] > 0.5) conv_partial[i] <- 1
}

fit_conv_partial <- quantile_spline(x, y, knots, tau = 0.5, convcons = conv_partial)
plot(x, y, pch = 16, cex = 0.5, col = "gray", main = "Partial Convexity")
lines(x_eval, spline_eval(fit_conv_partial, x_eval), col = "blue", lwd = 2)
abline(v = knots[6], col = "red", lty = 2)  # Boundary of constrained region

Mixed Constraints

Multiple shape constraints can be applied simultaneously.

# Increasing AND convex
fit_both <- quantile_spline(x, y, knots, tau = 0.5, 
                            monot = 1, convcons = 1)

plot(x, y, pch = 16, cex = 0.5, col = "gray", 
     main = "Mixed Constraints: Increasing + Convex")
lines(x_eval, spline_eval(fit_both, x_eval), col = "blue", lwd = 2)

Constraints with Multiple Knots

When using knot multiplicities, constraints need to be adapted to the reduced regularity. Constraints apply only on regular pieces : monotonicity is ensured on continuous pieces only.

# Create knots with multiplicity at a point
sn_mult <- c(0, 0, 0, 0, 0.3, 0.4,0.4,0.4, 0.5,0.5,0.5, 0.7, 1, 1, 1, 1)
#basis_mult <- Bspline_base(sn_mult, degree = 3)
knots_mult <- c(0, 0.3, 0.3,0.3,0.3,0.4, 0.5,0.7,0.7,0.7,0.7, 1)

# Fit with constraints
fit_mult <- quantile_spline(x, y, knots_mult, tau = 0.5, monot = 1)
plot(x, y, pch = 16, cex = 0.5, col = "gray", 
     main = "Constraints with Multiple knots")
lines(x_eval, spline_eval(fit_mult, x_eval), col = "blue", lwd = 2)

Recap: Visualizing Constraint Effects

par(mfrow = c(2, 2))

# No constraints
plot(x, y, pch = 16, cex = 0.4, col = "gray", main = "Unconstrained")
lines(x_eval, spline_eval(fit_uncon, x_eval), col = "red", lwd = 2)

# Monotonicity
plot(x, y, pch = 16, cex = 0.4, col = "gray", main = "Monotonicity")
lines(x_eval, spline_eval(fit_inc, x_eval), col = "blue", lwd = 2)

# Convexity
plot(x, y, pch = 16, cex = 0.4, col = "gray", main = "Convexity")
lines(x_eval, spline_eval(fit_convex, x_eval), col = "blue", lwd = 2)

# Mixed
plot(x, y, pch = 16, cex = 0.4, col = "gray", main = "Mixed")
lines(x_eval, spline_eval(fit_both, x_eval), col = "blue", lwd = 2)

par(mfrow = c(1, 1))

Karlin-Studden SOCP Formulation

The shape constraints are implemented using the Karlin-Studden characterization of non-negative polynomials of degree 2 or 3.

Depending on the degree of the spline (3 or 4), these constraints apply to the first or second derivative:

Constraints by Degree

Degree	Monotonicity	Convexity	Third Derivative
0 (constant)	+(constant)	X	X
1 (Linear)	+(linear)	X	X
2 (Quadratic)	+ (Linear)	+(constant)	X
3 (Cubic)	+(Karlin)	+(Linear)	+(linear (constant))
4 (Quartic)	+(Karlin)	+(Karlin)	+(linear)
5 (quintic)	X	+(Karlin)	Not implemented
\(\geq6\)	X	X	X

‘Karlin’ means in fact ‘quadratic’ inequality.

Other constraints pointwise constraints, at knots or not, are not implemented, but they would be implemented as simple linear inequalities.

The complete set of constraints forms a Second-Order Cone Program (SOCP). This SOCP formulation ensures exact shape constraints, unlike approximations that only enforce constraints at knots or at a finite set of points.

The particular structure of the problem allows the use of the CVXR package, which handles this type of mixed quadratic/linear conic problems efficiently through its DCP (Disciplined Convex Programming) framework and supports multiple solvers (CLARABEL, OSQP, ECOS, SCS).

Summary

The package provides: - Uniform constraints: Apply globally - Partial constraints: Apply on specific intervals - Mixed constraints: Multiple constraints simultaneously - SOCP formulation: Exact (not approximated) constraints on the whole intervals. Constraints concern monotone or convex