Optimal knot placement for fitting piecewise-continuous linear functions to a nonlinear function

Question

I encountered this problem in my research and it is turning out to be a surprisingly difficult one(for me, at least).

Suppose we have a univariate nonlinear function $f(x)$ where $x \in [L,U]$. Our goal is to approximate this nonlinear function with $n$ piecewise-continuous linear functions $g_{i}(x)$ within the given domain. We assume that $n$ is a pre-specified number. We define each line segment as follows: $$ g_{i}(x) = \frac{f(a_{i}) - f(a_{i-1})}{a_{i} - a_{i-1}} (x - a_{i-1}) + f(a_{i-1})\text{ for }a_{i-1} \leq x \leq a_{i} $$ where $a_{i}$ are knot points in $[L,U]$ and $i = 1,\ldots,n$. The first and the last knot points are fixed at the boundaries, that is, $a_{0} = L, a_{n} = U$. Also, the knot points are ordered and unique: $ a_{i} > a_{i-1}$ for $i=1,\ldots,n$.

I want to find the optimal placements for the knot points $a_{1},\ldots,a_{n-1}$, such that the overall squared-approximation error $e$ is minimized. We can pose the objective as follows: $$ \min_{a_{1},\ldots,a_{n-1}} \left\{ e = \int_{L}^{U} [f(x) - g_{i}(x)]^2 dx \right\} $$

This picture illustrates the problem: Piecewise Linear functions http://dl.dropbox.com/u/6809582/linearfunctions.png

The final optimization problem looks like the following (after a simple reformulation into a optimal-control-like form): $$ \begin{align*} &\min_{a_{1},\ldots,a_{n-1}} e(U)\\ s.t.\quad & \frac{de(x)}{dx} = [f(x) - g_{i}(x)]^2, \quad e(L) = 0\\ &g_{i}(x) = \frac{f(a_{i}) - f(a_{i-1})}{a_{i} - a_{i-1}} (x - a_{i-1}) + f(a_{i-1})\text{ for }a_{i-1} \leq x \leq a_{i}\\ & a_{0} = L, a_{n} = U\\ & a_{i} \geq a_{i-1} + \epsilon,\quad i=1,\ldots,n \end{align*} $$ This optimization problem is extremely difficult to solve numerically, owing to its nonsmoothness and nonconvexity.

Question: How do I solve this problem to global optimality? Can anyone provide any attacks (even partial ones)? Any simplifying properties?

Answer 1

In the limit of fine subdivision, the local goodness of fit depends on only two things: the absolute value of the local second derivative, and the local density of knots. For a segment of constant second derivative $a$ over an interval of length $c$, the integrated squared error over the interval comes out to $a^2c^5/120$. Given a small approximation interval (small enough that the second derivative does not vary by much over its length) whose squared error contribution is $E$, replacing that interval by two approximation intervals of half the width reduces the total squared error by approximately $E- 2(E/32)=(15/16)E$. The best segment to subdivide (according to this approximation) is thus whichever one contributes the most to the error. It follows that in the limit of many knots the optimum will have an equal contribution to the error coming from each segment. This is true, in the limit, when the local density of knots at $f(x)$ is proportional to $f^{\prime\prime}(x)^{2/5}$. A very good approximation to your optimal assignment problem can thus be obtained by making a plot, from $L$ to $U$, of $\int_L^x f^{\prime\prime}(s)^{2/5}\ ds$, dividing it into equally spaced horizontal strips, and placing your knots at the horizontal values where the integral plot crosses from one horizontal strip to the next.

Incidentally, if you keep the same number of straight-line approximations but allow them to cross the function (and back) in their interiors rather than at their endpoints, you can improve the integrated squared error by (in the limit) a factor of 6.

Answer 2

Do you really want to solve this problem, or do you just want to get a good approximation to the function? I would bet that your overall goal is actually to get a fast and reasonably accurate approximation to the function, and there are lots of ways of doing that without solving what you already know is a hard optimization problem.

If you have some good reason why you really need to solve the nonconvex nonsmooth optimization problem as you've formulated it, then consider using branch and bound with convex underestimating functions. If n is reasonably small, then a branch and bound code might actually be able to find a provable epsilon optimal solution in reasonable time.