• https://me.yahoo.com

# Univariate

Univeriate curve-fitting, interpolation, polynomial, spline, Akima

## Overview

Univariate interpolation is an area of curve-fitting which, as opposed to univariate regression analysis, finds the curve that provides an exact fit to a series of two-dimensional data points. It is called univariate as the data points are supposed to be sampled from a one-variable function. Compare this to multivariate interpolation, which aims at fitting data points sampled from a function of several variables.

Formally speaking, consider a series of $\inline&space;N$ data points $\inline&space;(x_1,&space;y_1),&space;(x_2,&space;y_2),&space;\ldots,&space;(x_N,&space;y_N)$ and, for the sake of simplicity, consider that $\inline&space;x_1&space;<&space;x_2&space;<&space;\ldots&space;<&space;x_N$, i.e. the points are distinct and are in increasing order with respect to $\inline&space;x$. By interpolating these data points we mean finding a function $\inline&space;f&space;:&space;[x_1,&space;x_n]&space;\to&space;\mathbb{R}$ such that:

$f(x_1)&space;=&space;y_1,&space;\quad&space;f(x_2)&space;=&space;y_2,&space;\quad&space;\ldots&space;\quad&space;f(x_N)&space;=&space;y_N.$

To make things even clearer, consider the following example. Let the data points be $\inline&space;(0,&space;0),&space;(\pi,&space;0),&space;(2&space;\pi,&space;0)$, shown in red in the image below.

As you may notice, we can choose $\inline&space;f(x)&space;=&space;\sin(x)$ to be the interpolation function, which is also shown in blue in the previous graph. It clearly satisfies the constraints of equation (1), but it's not necessarily the most straightforward choice. We might as well choose the constant function $\inline&space;f(x)&space;=&space;0$, which also satisfies the constraints of equation (1) and is therefore a valid interpolation function for the given data points.

As a conclusion, we are allowed to choose any type of function as long as it provides an exact fit to the given data points. Our choice should be based on some prior knowledge of the phenomena which generated that series of points.

## Polynomial Interpolation

In the following, let us assume that the interpolation function is polynomial, i.e. $\inline&space;f$ is of the type

$f(x)&space;=&space;a_0&space;+&space;a_1&space;x&space;+&space;a_2&space;x^2&space;+&space;\ldots&space;+&space;a_k&space;x^k,$

where $\inline&space;k$ is called the degree of $\inline&space;f$ and $\inline&space;a_0,&space;a_1,&space;\ldots,&space;a_k$ are some real numbers, called the coefficients of $\inline&space;f$. In order to find the expression for $\inline&space;f$, it suffices to find its coefficients. We are able to find them by writing the constraints of equation (1) in our particular case:

$\left\{\begin{array}{ccccccccccccc}&space;f(x_1)&space;&=&&space;y_1&space;&=&space;&a_0&space;&+&space;&a_1&space;x_1&space;&+&space;&a_2&space;x_1^2&space;&+&space;&\ldots&space;&+&space;&a_k&space;x_1^k&space;\\&space;f(x_2)&space;&=&&space;y_2&space;&=&space;&a_0&space;&+&space;&a_1&space;x_2&space;&+&space;&a_2&space;x_2^2&space;&+&space;&\ldots&space;&+&space;&a_k&space;x_2^k&space;\\&space;\ldots&space;&&&space;\ldots&space;&&&space;\ldots&space;\\&space;f(x_N)&space;&=&&space;y_N&space;&=&space;&a_0&space;&+&space;&a_1&space;x_N&space;&+&space;&a_2&space;x_N^2&space;&+&space;&\ldots&space;&+&space;&a_k&space;x_N^k&space;\\&space;\end{array}\right.$

This is a system of linear equations with $\inline&space;k&space;+&space;1$ unknowns, $\inline&space;a_0,&space;a_1,&space;\ldots,&space;a_k$ and $\inline&space;N$ equations. In order to have a unique solution, we need to require that $\inline&space;k&space;+&space;1&space;=&space;N$, or $\inline&space;k&space;=&space;N&space;-&space;1$; that is, we require that the degree of $\inline&space;f$ should be equal to the number of data points minus one. However, this is a worst case scenario since, by solving the above linear system, one may obtain $\inline&space;a_k&space;=&space;0$ which will decrease the degree of $\inline&space;f$ by one, and so on. As a general rule, the degree of $\inline&space;f$ will be strictly less than the number of given data points.

The function $\inline&space;f$ determined by solving the above system of linear equations is called the interpolation polynomial for the series of data points $\inline&space;(x_1,&space;y_1),&space;(x_2,&space;y_2),&space;\ldots,&space;(x_N,&space;y_N)$.

Directly solving this system of linear equations comes down to finding the inverse of its corresponding matrix, which might become computationally expensive. Another easier way of finding the coefficients of $\inline&space;f$ is by writing it in the Lagrange form

$f(x)&space;=&space;y_1&space;L_1(x)&space;+&space;y_2&space;L_2(x)&space;+&space;\ldots&space;+&space;y_N&space;L_N(x),$

where $\inline&space;L_i$ is defined through

$L_i(x)&space;=&space;&space;&space;\mathop{\prod_{j&space;=&space;1}^n}_{j&space;\neq&space;i}&space;&space;&space;\left(\frac{x&space;-&space;x_j}{x_i&space;-&space;x_j}\right),$

for all $\inline&space;i&space;=&space;1,&space;2,&space;\ldots,&space;N$. It can easily be shown that $\inline&space;L_i(x_i)&space;=&space;1$ and $\inline&space;L_i(x_j)&space;=&space;0$, for all $\inline&space;j&space;\neq&space;i$; therefore $\inline&space;f(x_i)&space;=&space;y_i$, so $\inline&space;f$ satisfies the constraints of equation (1). This provides a method of computing the interpolation polynomial $\inline&space;f$ known as Lagrange interpolation, which is also implemented as the component Interpolation/Lagrange.

Consider the function

$\psi(x)&space;=&space;\frac{\sin(2x)}{x}$

with $\inline&space;x$ in the interval $\inline&space;[\pi,&space;4&space;\pi]$ and let us sample 12 equidistant points from this function. Applying Lagrange interpolation on this series of data points results in the following approximation (shown in red):

## Runge's Phenomenon

When trying to estimate the error between the original function, from which the series of data points has been sampled, and the polynomial interpolation function $\inline&space;f$, one may notice the following phenomenon. Conside the Runge function:

$R(x)&space;=&space;\frac{1}{1&space;+&space;25&space;x^2}$

Now, consider a series of $\inline&space;N$ equidistant points $\inline&space;(x_1,&space;y_1),&space;(x_2,&space;y_2),&space;\ldots,&space;(x_N,&space;y_N)$ between $\inline&space;-1$ and $\inline&space;1$, where

$x_i&space;=&space;-1&space;+&space;(i&space;-&space;1)\frac{2}{N&space;-&space;1},\quad&space;y_i&space;=&space;f(x_i),$

for all $\inline&space;i&space;=&space;1,&space;2,&space;\ldots,&space;N$. Runge proved that the polynomial interpolation function corresponding to this set of data oscillates toward the end points of the interval $\inline&space;[-1,&space;1]$. More than that, the interpolation error at the ends of the interval tends to increase to infinity as you increase the number of equidistant data points. This is also known as Runge's phenomenon.

In the image below, the Runge function is depicted in red, while the 5-th degree (6 points) and 9-th degree (10 points) interpolation polynomials are shown in blue and green, correspondingly. Notice how the approximation error at the ends of the interpolation interval increases with the degree of $\inline&space;f$, or with the number of given equidistant points.

As a conclusion, polynomial interpolation in the case of a large number of equidistant points might generate large approximation errors toward the end points of the interpolation interval. This phenomenon can be avoided using spline interpolation, which is the subject of the next section.

## Spline Interpolation

Spline interpolation is somehow a generalization of polynomial interpolation, in that we do not necessarily have to find a single polynomial function to fit the data over the entire interval, but we rather try to find several polynomial functions to fit the data over each subinterval determined by two consecutive data points, while obeying some smoothness conditions. One of the advantages of this generalization is that the resulting interpolation function is less wiggly, as in the case of e.g. Lagrange interpolation.

Formally, consider the series of $\inline&space;N$ data points in the above paragraphs, which are distinct and ordered increasingly with respect to $\inline&space;x$. A spline interpolation function of degree $\inline&space;d&space;\geq&space;1$ for the given data points is a function $\inline&space;S&space;:&space;[x_1,&space;x_N]&space;\to&space;\mathbb{R}$ which satisfies the following conditions

• $\inline&space;S(x)&space;=&space;p_i(x)$, for all $\inline&space;x&space;\in&space;[x_i,&space;x_{i+1}]$ and all $\inline&space;i&space;=&space;1,&space;2,&space;\ldots,&space;N&space;-&space;1$, where $\inline&space;p_i$ is some polynomial function with degree less than or equal to $\inline&space;d$
• the derivatives of $\inline&space;S$ up to the order $\inline&space;n-1$ are all continuous in the given data points, which basically means that for all $\inline&space;i&space;=&space;1,&space;2,&space;\ldots,&space;N&space;-&space;1$ we require that:

$\left\{&space;\begin{array}{lcl}&space;p_i(x_{i+1})&space;&=&&space;p_{i+1}(x_{i+1})&space;\\&space;p_i'(x_{i+1})&space;&=&&space;p_{i+1}'(x_{i+1})&space;\\&space;p_i''(x_{i+1})&space;&=&&space;p_{i+1}''(x_{i+1})&space;\\&space;\ldots&&&space;\ldots&space;\\&space;p_i^{(n-1)}(x_{i+1})&space;&=&&space;p_{i+1}^{(n-1)}(x_{i+1}).&space;\end{array}&space;\right.$

Without getting into further technical details, we mention that the above relations lead to an under-determined system of linear equations. In order to obtain solution uniqueness to this system and completely determine the spline interpolation function $\inline&space;S$, $\inline&space;n&space;-&space;1$ more relations need to be set.

In the case of cubic splines ($\inline&space;n&space;=&space;3$), notice that we need $\inline&space;n&space;-&space;1&space;=&space;2$ additional relations to determine $\inline&space;S$. These are given, for example, by $\inline&space;p_1''(x_1)&space;=&space;p_2''(x_3)&space;=&space;0$. In this case, $\inline&space;S$ is called a natural cubic spline. Other kinds of cubic splines can be obtained using other pairs of conditions.

In the case when $\inline&space;n&space;=&space;1$, $\inline&space;S$ is piecewise linear and the method is called linear interpolation, which is also implemented as Interpolation/Linear. In the graph below you can see the way linear interpolation (in red) works in the case of 12 equidistant points sampled from the function $\inline&space;\psi$, previously defined.

In the case when $\inline&space;n&space;=&space;3$, $\inline&space;S$ is comprised of cubic polynomial functions on each subinterval and the method is called cubic interpolation, which is also implemented as Interpolation/Cubic. The following graph shows the cubic spline (in red) for the $\inline&space;\psi$ function.

Notice that higher degree spline interpolation results in smaller approximation errors. However, in order to escape from Runge's phenomenon, instead of increasing the degree of the spline, it is prefered to increase the number of data points. This can be achieved, for instance, by doing further experiments and obtaining further input data.

Akima interpolation is a particular type of third-degree spline interpolation, also implemented as Interpolation/Akima. Its main advantage is that, as opposed to cubic spline interpolation, it is applicable on successive intervals, so it does not require solving large systems of linear equations. Apart from being more computationally efficient, it also provides a more natural interpolation curve, closer to human intuition. The graph below shows the Akima interpolation curve for the $\inline&space;\psi$ function.

Lucian Bentea (August 2008)