# Polynomial regression interpolation? [duplicate]

Possible Duplicate:
Writing a function $f$ when $x$ and $f(x)$ are known

I'm not versed in mathematics, so you'll have to speak slowly...

If I want to fit a curve to the points,

X  Y
1  0.5
2  5.0
3  0.5
4  2.5
5  5.0
6  0.5


Where would I begin? For my purposes, this needs to be a sixth-order fit...

-

## marked as duplicate by robjohn♦, J. M., Phira, Zev ChonolesJun 11 '12 at 4:27

You're exactly right, Robjohn...That solved my problem. – Kevin Brown May 16 '12 at 14:03

First of all, since you have 6 points, you can only get a 5th order polynomial to fit. You need two points to define a line, 3 to define a quadratic, 4 for a cubic... n+1 for an nth-order.

The naive solution is simply to set up:

Y = A5 x^5 + A4 x^4 + A3 x^3 + A2 x^2 + A1 x + A0

For each data pair, plug in x, and generate an equation for the coefficients. You'll now have six equations in six unknowns, which is solvable, but computationally annoying if you're doing it by hand.

Example: for the second point, x=2, y=5

5 = A5 (2)^5 + A4 (2)^4 + A3 (2)^3 + A2 (2)^2 + A1 (2) + A0

and you end up with:

5 = 32 A5 + 16 A4 + 8 A3 + 4 A2 + 2 A1 + A0

Repeat for each data point, and then solve the system of equations for A0, A1, ..., A5

-
Good first step, Chris. How do I generate these equations? – Kevin Brown May 15 '12 at 20:17

Hopefully this isn't too far above your head. This is the most elegant way I know of to solve this, and it generalizes immediately to fitting of arbitrary degree and (with some knowledge of pseudoinverses) to best fit polynomials for overdetermined systems (as well as a fairly wide class of non-polynomial fits).

Take the Vandermonde matrix: $A= \left( \begin{array}{ccc} 1 & 1 & 1 & 1 & 1 & 1 \\ 1 & 2 & 2^2 & 2^3 & 2^4 & 2^5 \\ 1 & 3 & 3^2 & 3^3 & 3^4 & 3^5 \\ 1 & 4 & 4^2 & 4^3 & 4^4 & 4^5 \\ 1 & 5 & 5^2 & 5^3 & 5^4 & 5^5 \\ 1 & 6 & 6^2 & 6^3 & 6^4 & 6^5 \\ \end{array} \right)$

and the column vector $b=\left( \begin{array}{ccc} f(1) & f(2) & f(3) & f(4) & f(5) & f(6)\\ \end{array} \right)^{T}$, where $T$ denotes the transpose and $f$ is the undetermined polynomial of degree 5 (or less).

The equation you want to solve is the linear equation $A y = b$ for a column vector $y$. The reason for this is that if you take an arbitrary row of the matrix and multiply it by $y$, what you get is $f(n) = a_0 +a_1 n + a_2 n^2 + a_3 n^3 + a_4 n^4 + a_5 n^5$, where $y=\left( \begin{array}{ccc} a_0 & a_1 & a_2 & a_3 & a_4 & a_5\\ \end{array} \right)^{T}$ and $n \in \{1,2,3,4,5,6\}$.

The solution is, of course, $y = A^{-1}b$. The general formula for the inverse of a Vandermonde matrix is known, and in this case the inverse is:

$A^{-1} = \frac{1}{120} \left( \begin{array}{ccc} 720 & -1800 & 2400 & -1800 & 720 & -120 \\ -1044 & 3510 & -5080 & 3960 & -1620 & 274 \\ 580 & -2305 & 3720 & -3070 & 1300 & -225 \\ -155 & 685 & -1210 & 1070 & -475 & 85 \\ 20 & -95 & 180 & -170 & 80 & -15 \\ -1 & 5 & -10 & 10 & -5 & 1 \\ \end{array} \right)$

(Please check this if you intend to use it as I could have easily made a typo)

Now with simple multiplication you can find $a_0, \ldots, a_5$, and the polynomial you want is $f(x) = a_0 + a_1 x + a_2 x^2 + \cdots + a_5 x^5$.

EDIT I feel it necessary to point out that this isn't a terribly good method for interpolation in most real-world situations, and it's absolutely terrible for extrapolation. See, for example, Runge's Phenomenon. If you are going to use this method for interpolation, you should use the lowest degree polynomials that accurately represent the data. There are methods that don't have the same issues (though of course no method is perfect); you should ask a separate question about these if you are interested.

-
Dear Logan, it is even better to remember the LDU-decomposition of the vandermonde matrix intro triangular and diagonal factors and their triangular inverses; that have then entries independent of the dimension and are constant for any dimension.(We get the Pascal-matrix and the matrix of Stirling numbers of 2nd resp 1st kind): $\small V=L\cdot D \cdot U = P \cdot F \cdot S2~$ and $\small V^{-1}=U^{-1}\cdot D^{-1} \cdot L^{-1} = S1~ \cdot F^{-1} \cdot P^{-1}$ for any dimension. – Gottfried Helms May 16 '12 at 1:15
This is all good theory. Thanks! How can I do this in excel or C#...? – Kevin Brown May 16 '12 at 18:33
Both Excel and C# have routines for matrix multiplication and inversion. I'm not familiar with either, though. In excel, to find $A^{-1}$ use the MINVERSE function: addictivetips.com/microsoft-office/…. To do multiplication, use MMULT: addictivetips.com/microsoft-office/…. You'll need to enter $b$ by hand. You can probably write a VBA script so that you don't have to enter all the entries of $A$. Since this isn't a programming board, maybe ask at StackOverflow if you still have trouble. – Logan Maingi May 16 '12 at 18:43

Here is a step by step way of finding a polynomial $f(x)$ such that for $1 \leq i \leq 6$, $f(i)$ has the values you stated. We will "build up" $f(x)$ slowly, calling the intermediate stages $f^{(1)}(x), f^{(2)}(x), \ldots$ and so on.

• Set $y = f^{(1)}(x) = 0.5$ where the right side doesn't really depend on $x$ and the $0.5$ was chosen from the given data: when we set $x=1$ in $f^{(1)}(x)$, this (zero-degree) polynomial evaluates to $0.5$, and thus fits the point $(1,0.5)$ that was given to you.

• If we set $x=2$ in $f^{(1)}(x)$, we will get $0.5$ but we really need to get $5.0$. So, let us test $$y = f^{(2)}(x) = f^{(1)}(x) + c(x-1) = 0.5 + c(x-1)$$ where $c$ is just some number (for now) to see what happens. If we set $x=1$ in $f^{(2)}(x)$, we get $0.5$ as before (yay!) but if we set $x=2$, we get $0.5 + c$ and so by choosing $c = 4.5$, we can get $f^{(2)}(x)$ to give the desired values for both $x = 1$ and $x = 2$.

$f^{(2)}(x) = 0.5 + 4.5(x-1)$ is a polynomial such that $f^{(2)}(1) = 0.5$ and $f^{(2)}(2) = 5.0$.\ Notice that $4.5 = 5.0 - f^{(1)}(2)$ which I will write for future reference as $\displaystyle\frac{5.0 - f^{(1)}(2)}{2-1}$.

• $f^{(2)}(3) = 0.5 + 4.5(3-1) = 9.5$ which is not what we want. So, let us try $$y = f^{(3)}(x) = f^{(2)}(x) + c(x-1)(x-2)$$ where, as before, we will choose the value of $c$ in just a bit. Note that if we set $x$ to either $1$ or $2$ in $f^{(3)}(x)$, that last term evaluates to $0$ and so $f^{(3)}(1) = f^{(2)}(1) = 0.5$, $f^{(3)}(2) = f^{(2)}(2) = 5.0$ and so we are OK there. But, $$f^{(3)}(3) = f^{(2)}(3) + c(2)(1) = 9.5 + 2c$$ and so choosing $c = -4.5$ will make $f^{(3)}(3)$ evaluate to the desired $0.5$.

$f^{(3)}(x) = 0.5 + 4.5(x-1) - 4.5(x-1)(x-2)$ is a polynomial such that $f^{(3)}(1) = 0.5$, $f^{(3)}(2) = 5.0$, and $f^{(3)}(3) = 0.5$. Note that $-4.5 = \displaystyle\frac{0.5 - f^{(2)}(3)}{(3-1)(3-2)}$.

Since I am getting tired of typing, let me briefly explain how the rest of the thing will work. We set $$f^{(4)}(x) = f^{(3)}(x) + \frac{2.5 - f^{(3)}(4)}{(4-1)(4-2)(4-3)}(x-1)(x-2)(x-3)$$ where the second term vanishes when we set $x = 1, 2, 3$ and so we get the desired values of $y$, and when we set $x = 4$, we get $$f^{(4)}(4) = f^{(3)}(4) + \frac{2.5 - f^{(3)}(4)}{(4-1)(4-2)(4-3)}(4-1)(4-2)(4-3) = 2.5$$ which is exactly what we want. Lather, rinse, and repeat for $f^{(5)}(x)$ and $f^{(6)}(x)$.