Introduction
Sometimes you don't have data where you need it. Maybe you have figures collected sporadically but you need to present them as if they were collected at regular intervals, say once a week. Or maybe you need to make your best guess at missing data points. This is a job for interpolation.
Sometimes instead of wanting to guess at output values for inputs that you choose, you want to do the opposite: You have an output in mind and you need to guess what input would give you the output that you want. This calls for inverse interpolation.
This article presents the simplest and most widely applicable forms of interpolation — linear interpolation and quadratic interpolation. The article concludes with some recommendations for what to try next when these forms of interpolation are not good enough.
Linear Interpolation
There are many ways to fill in missing data, ranging from very simple to very sophisticated. Linear interpolation is the simplest method. It is also one of the most robust methods, i.e. it is likely to give reasonable answers under a wide variety of circumstances, including when substantial noise is added to the inputs. Unfortunately linear interpolation is also the least accurate method. However, accuracy is overrated. Or at least robustness is underrated. Robustness is often more important than extra accuracy, especially if you're not exactly confident that you know what you're doing. I would give the same advice for interpolation that the Extreme Programming folks give for software development: first try the simplest thing that could possibly work.
If you have two inputs, x_{0}
and x_{1}
, and two corresponding outputs y_{0}
and y_{1}
, the equation of the line connecting (x_{0}, y_{0})
and (x_{1}, y_{1})
is the following:
y = y<sub>0</sub>(x  x<sub>1</sub>)/(x<sub>0</sub>  x<sub>1</sub>) + y<sub>1</sub>(x  x<sub>0</sub>)/(x<sub>1</sub>  x<sub>0</sub>)
With this formula, you can stick in any value of x
you want and get out a new value of y
. So if you have a value x_{2}
and you want to guess its corresponding output y_{2}
, then you have this equation:
y<sub>2</sub> = y<sub>0</sub>(x<sub>2</sub>  x<sub>1</sub>)/(x<sub>0</sub>  x<sub>1</sub>) + y<sub>1</sub>(x<sub>2</sub>  x<sub>0</sub>)/(x<sub>1</sub>  x<sub>0</sub>)
Linear interpolation is most reliable if the x
you stick in is between the values x_{0}
and x_{1}
. If x
is not between these two values, technically you are extrapolating rather than interpolating. This still works as long as the new value of x
isn't too far from x_{0}
and x_{1}
. The further the new x
value is from the input values used to specify the line, the more suspicious you should be of your output.
Now suppose you have the data points (x_{0}, y_{0})
and (x_{1}, y_{1})
, but instead of trying to predict a new y
value you want to predict a new x
value. That is, you have a y
value you're trying to get out and you want to guess what input x
would give you that output. Then you can reverse the roles of x
and y
in the equation above and get the following:
x<sub>2</sub> = x<sub>0</sub>(y<sub>2</sub>  y<sub>1</sub>)/(y<sub>0</sub>  y<sub>1</sub>) + x<sub>1</sub>(y<sub>2</sub>  y<sub>0</sub>)/(y<sub>1</sub>  y<sub>0</sub>)
As before, this works best if the new value y_{2}
is between the previous y
values y_{0}
and y_{1}
. If it's not between these values but close to one of them, the result is likely to be useful. The further out y
gets, the more suspicious you should be of the result.
The formulas above are mathematically correct for any values you stick in. When I say that you should be suspicious of the output in extreme circumstances, it's not because some approximation is going on. However, in practice, the assumption that the points (x_{0}, y_{0})
and (x_{1}, y_{1})
can be used to accurately predict the point (x_{2}, y_{2})
may not hold when x_{2}
or y_{2}
are far from their predecessors.
When you write code to do linear interpolation, the only thing to be careful about is input validation. For (ordinary) interpolation, it's important to verify that x_{1}
does not equal x_{2}
so that the interpolation function doesn't divide by zero. Given arrays of x
and y
values, the following code fits a straight line to (x[0], y[0])
and (x[1], y[1])
and uses this line to predict a y[2]
value for x[2]
.
if (x[0] == x[1]  x[0] == x[2]  x[1] == x[2])
// report error
else
y[2] = y[0]*(x[2]  x[1])/(x[0]  x[1]) + y[1]*(x[2]  x[0])/(x[1]  x[0]);
Similarly, for inverse interpolation you need to make sure your y
values are distinct so that you don't divide by zero. The following uses a straight fitted line to (x[0], y[0])
and (x[1], y[1])
to predict the x[2]
value for y[2]
.
if (y[0] == y[1]  y[0] == y[2]  y[1] == y[2])
// report error
else
x[2] = x[0]*(y[2]  y[1])/(y[0]  y[1]) + x[1]*(y[2]  y[0])/(y[1]  y[0]);
If you just need to do a quick interpolation on a small amount of data rather than write a program that does interpolation, you may want to use an online linear interpolation calculator. This interpolator is implemented in handwritten clientside JavaScript and so you can read the source.
Quadratic Interpolation
As the Extreme Programming folks would recommend, when the simplest thing doesn't work, try the next simplest thing that could possibly work. For interpolation, that means quadratic interpolation. Instead of fitting a straight line to two points, quadratic interpolation fits a parabola to three points. To figure out how to generalize the formulas above to quadratics, look back at the equation for linear interpolation. The term (x  x_{1})/(x_{0}  x_{1})
is 1
when x = x_{0}
and 0
when x = x_{0}
. Therefore the term y_{0}(x  x_{1})/(x_{0}  x_{1})
is y_{0}
at x_{0}
and 0
at x_{1}
. Similarly, (x  x_{0})/(x_{1}  x_{0})
is 1
when x = x_{1}
and 0
when x = x_{0}
, and so y_{1}(x  x_{0})/(x_{1}  x_{0})
is y_{1}
at x_{1}
and 0
at x_{0}
.
For quadratic interpolation, we follow a similar pattern, constructing quadratic polynomials that are 1
at one of the given x
s and 0
at the others. Then when we multiply by the corresponding y
values and add up terms, the result has the necessary y
values for each x
. So the quadratic polynomial fitting the points (x_{0}, y_{0})
, (x_{1}, y_{0})
, and (x_{2}, y_{2})
is:
y<sub>0</sub>P<sub>0</sub>(x) + y<sub>1</sub>P<sub>1</sub>(x) + y<sub>2</sub>P<sub>0</sub>(x)
where:

P<sub>0</sub>(x) = (x  x<sub>1</sub>)(x  x<sub>2</sub>)/((x<sub>0</sub>  x<sub>1</sub>)(x<sub>0</sub>  x<sub>2</sub>))

P<sub>1</sub>(x) = (x  x<sub>0</sub>)(x  x<sub>2</sub>)/((x<sub>1</sub>  x<sub>0</sub>)(x<sub>1</sub>  x<sub>2</sub>))

P<sub>2</sub>(x) = (x  x<sub>0</sub>)(x  x<sub>1</sub>)/((x<sub>2</sub>  x<sub>0</sub>)(x<sub>2</sub>  x<sub>1</sub>))
The polynomials P_{i}
above are called "Lagrange" polynomials. For cubic and higher interpolation, the pattern is the same: first construct polynomials that are 1
at one of the x
s and 0
at all other x
s, then multiply each by the corresponding y
values and sum.
If your inputs are free of noise, quadratic interpolation can give much better accuracy than linear interpolation. For example, in mathematical tables, the given values are precise to many decimal places, but you may be interested in a value not in the table. Say a function is tabulated at 0.1
, 0.2
and 0.3
but you want to know its value at 0.17
. You could probably get much better accuracy using a parabola to fit all three points rather than using a line to just fit the first two points. On the other hand, if your inputs have a substantial amount of error, some sort of random noise, then quadratic interpolation could magnify that noise by overreacting to the noise.
To implement the above equations in software, we need to verify that all the x
values are distinct in order to avoid dividing by zero.
if (x[0] == x[1]  x[0] == x[2]  x[1] == x[2]  x[0] == x[3] 
x[1] == x[3]  x[2] == x[3])
{
// report error;
}
else
{
y[3] = y[0]*(x[3]  x[1])*(x[3]  x[2])/((x[0]  x[1])*(x[0]  x[2]));
y[3] += y[1]*(x[3]  x[0])*(x[3]  x[2])/((x[1]  x[0])*(x[1]  x[2]));
y[3] += y[2]*(x[3]  x[0])*(x[3]  x[1])/((x[2]  x[0])*(x[2]  x[1]));
}
As with the linear interpolator, there is an online quadratic interpolation calculator implemented with clientside JavaScript.
Now what if we want to use interpolation to solve for an x
value? Say we have three points (x_{0}, y_{0})
, (x_{1}, y_{0})
, and (x_{2}, y_{2})
and we want to solve for the x_{3}
value corresponding to a given y_{3}
? We could fit a polynomial to the first three points exactly as above and solve for x
. That would generally not be a good idea. Instead, we reverse the roles of x
and y
and imagine y
as the independent variable.
Note that this is where quadratic interpolation differs from linear interpolation. With linear interpolation, reversing the roles of x
and y
is the same as fitting first as a function of x
and then solving for a missing x
. With quadratic interpolation, the analogous steps are not the same. Fitting the points as a function of x
first then solving the resulting equation would amount to finding the roots of a quadratic equation. This might not be possible, or it might give two different answers. Even if it gives a single answer, that answer might amplify errors in the data. On the other hand, by simply reversing the roles of x
and y
, we have a simple, wellbehaved solution.
What If Linear and Quadratic Interpolation Aren't Good Enough?
In theory, you could fit higher order polynomials. You could fit a third degree polynomial to four points, or a fourth degree polynomial to give points, etc. This is generally not a good idea. Fitting higher degree polynomials amplifies errors in the data.
If quadratic interpolation isn't good enough, you may need some more sophisticated form of interpolation. For example, natural cubic splines are useful in many contexts. However, it's hard to say much in general. Beyond simple linear or quadratic interpolation, the best technique depends heavily on the problem context. You may need to abandon interpolation entirely and use a DSP algorithm, such as a lowpass filter, or do some sort of statistical regression. A good place to start looking for more information would be the Numerical Recipes book.
History
 11^{th} September, 2008: Initial post
 8^{th} October, 2008: Rewritten to include quadratic interpolation and code samples