# Lagrange Interpolation

## What it is

Lagrange Interpolation allows for the reconstruction of a polynomial with the help of a set of values and their respective evaluations on such polynomial.

More specifically, a degree $n$ single-variable polynomial $f(x)$ can be uniquely reconstructed via $n+1$ polynomial evaluations, that is $n+1$ pairs of points of the form $(x, f(x))$.

## How it works

In this writeup we’ll solely focus on single-variable polynomials which are defined as follows:

$f(x) = a_0 + a_1x^1 + a_2x^2 + ... + a_nx^n$

Where $n$ is the polynomial’s degree and its coefficients are $a_0$, $a_1$, $a_2$, … $a_n$.

To evaluate such a polynomial at $x$, one substitutes the $x$-values and carries out the necessary math to derive what’s usually referred to as the $y$-value.

Given that we have $x$- and $y$ values, we can interpret them as coordinates and plot them in a 2-dimensional Euclidean Space using the Cartesian Coordinate System which allows for the visualization of the polynomial.

If you’ve studied polynomials in school, you might remember that a degree $1$ polynomial represents a line, a degree $2$ polynomial a parabola, a degree $3$ polynomial a cubic function, etc.

Given the definition of a polynomial, we can now state our goal as being able to reconstruct any degree $n$ polynomial by having access to at least $n+1$ pairs of $x$- and $y$ values (which we’ll call points from now on).

### System of Linear Equations

Before we dive straight into the mechanics of Lagrange Interpolation we should motivate its versatility and usefulness by looking at three different examples where we reconstruct a degree $0$, degree $1$ and degree $2$ polynomial without Lagrange Interpolation.

#### Example 1: Degree 0 Polynomial

Our first example is the “trivial” example in which we have a degree $0$ polynomial we want to reconstruct.

We therefore have $n=0$ and need $n + 1 = 0 + 1 = 1$ point.

In our case, the point we’re using for the evaluation is $(1, 2)$, so $f(1) = 2$.

Looking at the definition of a polynomial above, we can see that our degree $0$ polynomial only consists of a single coefficient $a_0$ which results in $f(x) = a_0$.

Given that we know that $f(1) = 2$ we can derive that $a_0 = 2$, so we’ve successfully reconstructed the polynomial as:

$f(x) = 2$

The following graph shows the polynomial being drawn in a two-dimensional coordinate system (note the point we’ve used to reconstruct the polynomial):

#### Example 2: Degree 1 Polynomial

To make things more interesting we’ll now be looking at a degree $1$ polynomial to reconstruct (a line).

In this case $n=1$ and we need $n + 1 = 1 + 1 = 2$ points for reconstruction.

Assuming our points are $(1, 2)$ and $(3, 4)$, this means that $f(1) = 2$ and $f(3) = 4$.

As per the definition of a polynomial, we’re looking for a degree $1$ polynomial $f(x) = a_0 + a_1x^1$.

To find the polynomial in question we can create a system of linear equations using the polynomial “skeleton” and the values of our points. Solving this linear system (e.g. with Gaussian Elimination) allows us to reconstruct our degree $1$ polynomial.

\begin{aligned} a_0 + a_1(1)^1 &= 2 \\ a_0 + a_1(3)^1 &= 4 \\ \\ a_0 + 1a_1 &= 2 && \mid -1a_1 \\ a_0 + 3a_1 &= 4 \\ \\ a_0 &= 2 - 1a_1 \\ \\ a_0 + 3a_1 &= 4 \\ 2 - 1a_1 + 3a_1 &= 4 \\ 2 + 2a_1 &= 4 && \mid -2 \\ 2a_1 &= 2 && \mid \div 2 \\ a_1 &= 1 \\ \\ a_0 + 1a_1 &= 2 \\ a_0 + 1 \times 1 &= 2 \\ a_0 + 1 &= 2 && \mid -1 \\ a_0 &= 1 \end{aligned}

Given that we isolated the $a_0$ and $a_1$ values, we can reconstruct the polynomial as:

\begin{aligned} f(x) &= a_0 + a_1x^1 \\ &= 1 + 1x \\ &= 1 + x \end{aligned}

You might remember a saying from your math teacher that goes something like this: “A line is uniquely defined via two points”. The math we did above serves as a proof of such statement.

Our degree $1$ polynomial can be drawn in a two-dimensional coordinate system like follows (note the two points we’ve used for reconstruction):

#### Example 3: Degree 2 Polynomial

So far we’ve successfully reconstructed polynomials of degree $0$ and $1$. As a last step we’ll now focus on the reconstruction of a degree $2$ polynomial (a parabola) which we’ll reuse throughout the rest of this article.

We have $n=2$ and therefore need $n + 1 = 2 + 1 = 3$ points to reconstruct the polynomial successfully.

Our points for this example are $(1, 2)$, $(2, 3)$ and $(3, 6)$ which translates into $f(1) = 2$, $f(2) = 3$ and $f(3) = 6$.

Checking the definition of a polynomial we can see that we’re looking for a degree $2$ polynomial of the form $f(x) = a_0 + a_1x^1 + a_2x^2$.

Once again, we can use our “skeleton” polynomial and our points to setup a system of linear equations. We then solve this system with Gaussian Elimination to reconstruct our degree $2$ polynomial.

\begin{aligned} a_0 + a_1(1)^1 + a_2(1)^2 &= 2 \\ a_0 + a_1(2)^1 + a_2(2)^2 &= 3 \\ a_0 + a_1(3)^1 + a_2(3)^2 &= 6 \\ \\ a_0 + 1a_1 + 1a_2 &= 2 \\ a_0 + 2a_1 + 4a_2 &= 3 && \mid \times 9 \\ a_0 + 3a_1 + 9a_2 &= 6 && \mid \times (-4) \\ \\ a_0 + 1a_1 + 1a_2 &= 2 && \mid \times (-4) \\ a_0 + 2a_1 + 4a_2 &= 3 && \mid \times 1 \\ 5a_0 + 6a_1 + \cancel{0a_2} &= 3 \\ \\ a_0 + 1a_1 + 1a_2 &= 2 \\ -3a_0 - 2a_1 + \cancel{0a_2} &= -5 && \mid \times 3 \\ 5a_0 + 6a_1 + \cancel{0a_2} &= 3 && \mid \times 1 \\ \\ a_0 + 1a_1 + 1a_2 &= 2 \\ -3a_0 - 2a_1 + \cancel{0a_2} &= -5 \\ -4a_0 + \cancel{0a_1} + \cancel{0a_2} &= -12 \\ \\ -4a_0 &= -12 && \mid \div (-4) \\ a_0 &= 3\\ \\ -3a_0 - 2a_1 &= -5 \\ -3(3) - 2a_1 &= -5 \\ -9 -2a_1 &= -5 && \mid +9 \\ -2a_1 &= 4 && \mid \div (-2) \\ a_1 &= -2 \\ \\ a_0 + 1a_1 + 1a_2 &= 2 \\ 3 + 1(-2) + 1a_2 &= 2 \\ 3 - 2 + 1a_2 &= 2 \\ 1 + 1a_2 &= 2 && \mid -1 \\ 1a_2 &= 1 \\ a_2 &= 1 \end{aligned}

With values for $a_0$, $a_1$ and $a_2$ we can reconstruct the polynomial as:

\begin{aligned} f(x) &= a_0 + a_1x^1 + a_2x^2 \\ &= 3 - 2x + 1x^2 \\ &= 3 - 2x + x^2 \end{aligned}

Drawing this degree $2$ polynomial in a two-dimensional coordinate system results in the following diagram (again, note the points we’ve used to reconstruct the polynomial):

### Lagrange Interpolation

As we’ve seen with our examples above, we can use techniques such as a system of linear equations to reconstruct a polynomial reliably. However, as the degree of the polynomial increases, so does the complexity of our linear systems.

There’s got to be a better way to reconstruct polynomials and luckily enough there is! We can use a tool called Lagrange Interpolation to reconstruct any degree $n$ polynomial with the help of $n+1$ points without reaching for a system of linear equations.

The main observation that mathematicians like Joseph-Louis Lagrange had is that any degree $n$ polynomial can be rewritten as a sum of $n+1$ “partial polynomials” which are usually referred to as “basis polynomials”. Such basis polynomials need to assume a non-zero $y$-value for their respective $x$-value and a zero $y$-value for all the other $n$ $x$-values.

#### Basis Polynomial Construction

To better understand what this definition means, let’s take a look at our degree $2$ polynomial we’ve worked with before. Its points are $(1, 2)$, $(2, 3)$ and $(3, 6)$.

So, for every basis polynomial $\mathcal{L_i}(x)$ that takes the $x$-value of one of the three points as input, we want it to behave in the following way:

$\begin{equation*} \mathcal{L_i}(x) = \begin{cases} 1 & x = x_i \\ 0 & x = x_j && j \ne i \end{cases} \end{equation*}$

Taking the point $(1, 2)$ as an example, we want the basis polynomial $\mathcal{L_1}(x)$ defined for this point to assume the following values:

$\mathcal{L_1}(1) = 2$ $\mathcal{L_1}(2) = 0$ $\mathcal{L_1}(3) = 0$

The basis polynomial $\mathcal{L_2}(x)$ we define for point $(2, 3)$ should assume the following values:

$\mathcal{L_2}(1) = 0$ $\mathcal{L_2}(2) = 3$ $\mathcal{L_2}(3) = 0$

And the basis polynomial $\mathcal{L_3}(x)$ that’s defined for the point $(3, 6)$ should output the following values:

$\mathcal{L_3}(1) = 0$ $\mathcal{L_3}(2) = 0$ $\mathcal{L_3}(3) = 6$

Thinking about this a bit, we can utilize the observation that an equation like $(x - a)(x - b)$ is $0$ if either $x = a$ or $x = b$, but non-zero otherwise. If we have a fraction we could use this equation in the numerator in conjunction with a proper denominator to either turn the whole fraction into a $0$ (a fraction with a $0$ in the numerator will always be $0$) or a $1$. We could then use this “binary switch” together with a rescaling by the respective $y$-value to compute the correct value.

Equipped with this knowledge we can define our basis polynomials as follows:

$\mathcal{L_1}(x) = \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)}$ $\mathcal{L_2}(x) = \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)}$ $\mathcal{L_3}(x) = \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)}$

To show that these basis polynomials work (almost) as desired, we can take a look at $\mathcal{L_1}(x)$ and evaluate it at the $x$-values $1$, $2$ and $3$ (the basis polynomials $\mathcal{L_2}(x)$ and $\mathcal{L_3}(x)$ behave similarly which we won’t outline here).

$\mathcal{L_1}(1) = \frac{(1 - 2)(1 - 3)}{(1 - 2)(1 - 3)} = \frac{2}{2} = 1$ $\mathcal{L_1}(2) = \frac{(2 - 2)(2 - 3)}{(1 - 2)(1 - 3)} = \frac{0}{2} = 0$ $\mathcal{L_1}(3) = \frac{(3 - 2)(3 - 3)}{(1 - 2)(1 - 3)} = \frac{0}{2} = 0$

#### Rescaling of Basis Polynomials

As we can see, we were able to produce a polynomial that assumes a value of either $0$ or $1$ but it doesn’t assume the correct $y$-value yet. This however is easy to fix as we can simply multiply the result of our basis polynomial $\mathcal{L_i}(x)$ by the respective $y$-value as we’ve alluded to above.

We define these rescaled basis polynomials as follows:

\begin{aligned} P_1(x) &= y_1 \times \mathcal{L_1}(x) \\ &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} \end{aligned} \begin{aligned} P_2(x) &= y_2 \times \mathcal{L_2}(x) \\ &= y_2 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} \end{aligned} \begin{aligned} P_3(x) &= y_3 \times \mathcal{L_3}(x) \\ &= y_3 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \end{aligned}

Doing a quick sanity check, we can confirm that this rescaling does the trick:

\begin{aligned} P_1(1) &= y_1 \times \mathcal{L_1}(x) \\ &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} \\ &= 2 \times \frac{(1 - 2)(1 - 3)}{(1 - 2)(1 - 3)} \\ &= 2 \times 1 \\ &= 2 \end{aligned} \begin{aligned} P_1(2) &= y_1 \times \mathcal{L_1}(x) \\ &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} \\ &= 2 \times \frac{(2 - 2)(2 - 3)}{(1 - 2)(1 - 3)} \\ &= 2 \times 0 \\ &= 0 \end{aligned}

#### Summation of Rescaled Basis Polynomials

Given that we’ve now defined these rescaled basis polynomials that exhibit the desired behavior of assuming the correct $y$-value given a point’s $x$-value and returning $0$ as a $y$-value for all the other $n$ point’s $x$-values, we can now add them together to create our final polynomial $P$ which is called the “interpolation polynomial”:

\begin{aligned} P(x) &= P_1(x) + P_2(x) + P_3(x) \\ &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + y_2 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + y_3 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \end{aligned}

To prove that this is indeed the polynomial that describes our parabola, we can substitute the $y$-values and apply some basic algebra:

\begin{aligned} P(x) &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + y_2 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + y_3 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= 2 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + 3 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + 6 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= 2 \times \frac{x^2 - 3x - 2x + 6}{1 - 3 - 2 + 6} + 3 \times \frac{x^2 - 3x - x + 3}{4 - 6 - 2 + 3} + 6 \times \frac{x^2 - 2x - x + 2}{9 - 6 - 3 + 2} \\ &= \frac{\cancel{2}(x^2 - 3x - 2x + 6)}{\cancel{2}} + \frac{3(x^2 - 3x - x + 3)}{-1} + \frac{\cancel{6}(x^2 - 2x - x + 2)}{\cancel{2}} \\ &= (x^2 - 5x + 6) - 3(x^2 - 4x + 3) + 3(x^2 - 3x + 2) \\ &= x^2 - 5x + 6 - 3x^2 + 12x - 9 + 3x^2 - 9x + 6 \\ &= x^2 - 2x + 3 \\ &= 3 - 2x + x^2 \end{aligned}

Comparing this result with the result from above where we’ve used a system of linear equations, we can indeed see that we reconstructed the correct polynomial! This time around however, we’ve summed up rescaled basis polynomials to get an interpolation polynomial. A process that’s called Lagrange Interpolation.

#### Generalization for Arbitrary Degree $n$ Polynomials

Now that we’ve discovered the mechanics behind Lagrange Interpolation from first principles and went through an instantiation of such, we can look at its building blocks once again to define a formula we can use to reconstruct arbitrary degree $n$ polynomials.

Starting with the basis polynomials $\mathcal{L_i}(x)$, we can generalize its construction as follows:

\begin{aligned} \mathcal{L_i}(x) &= \frac{x - x_j}{x_i - x_j} \\ &= \frac{x - x_j}{x_i - x_j} \times \frac{x - x_{j+1}}{x_i - x_{j + 1}} \times ... \times \frac{x - x_n}{x_i - x_n} && \text{for j \ne i} \\ &= \prod_{\substack{j=0 \\ j \ne i}}^{n} \frac{x - x_j}{x_i - x_j} \end{aligned}

Next up, we can define the rescaled basis polynomial $P_i(x)$ as:

$P_i(x) = y_i \times \mathcal{L_i}(x)$

And finally our interpolation polynomial $P(x)$ can be defined as the sum of all the rescaled basis polynomials:

\begin{aligned} P(x) &= P_i(x) + P_{i + 1}(x) + ... + P_n(x) \\ &= \sum_{i = 0}^{n} P_i(x) \\ &= \sum_{i = 0}^{n} y_i \times \mathcal{L_i}(x) \\ &= \sum_{i = 0}^{n} y_i \left(\prod_{\substack{j=0 \\ j \ne i}}^{n} \frac{x - x_j}{x_i - x_j}\right) \end{aligned}

And so the following formula is what’s known as- and used for Lagrange Interpolation:

$P(x) = \sum_{i = 0}^{n} y_i \left(\prod_{\substack{j=0 \\ j \ne i}}^{n} \frac{x - x_j}{x_i - x_j}\right)$

## Why it works

The key takeaway of the Lagrange Interpolation algorithm is that a complex degree $n$ polynomial can be uniquely represented by summing up $n+1$ “partial polynomials” (which are called basis polynomials) that each individually contribute to the polynomial which should be reconstructed.

This “contribution” of each individual basis polynomial is defined as given its corresponding point’s $x$-value it should assume the correct $y$-value, but should return the $y$-value of $0$ whenever the $x$-value of one of the remaining $n$ points is supplied.

Using $x$-values we can therefore “toggle on” and “toggle off” the desired $y$-values as well as calculate the correct $y$-values for $x$-values that lie between two points.

\begin{aligned} P(x_i) &= y_0\mathcal{L_0}(x_i) + y_1\mathcal{L_1}(x_i) + ... + y_i\mathcal{L_i}(x_i) + ... + y_n\mathcal{L_n}(x_i) \\ &= 0 + 0 + y_i + 0 \\ &= y_i \end{aligned}

Using our example from above once again, we can prove that we’re indeed able to “enable” / “disable” the correct $y$-values given their corresponding $x$-value:

\begin{aligned} P(x) &= P_1(x) + P_2(x) + P_3(x) \\ &= y_1 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + y_2 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + y_3 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= 2 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + 3 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + 6 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \end{aligned} \begin{aligned} P(1) &= 2 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + 3 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + 6 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= (2 \times 1) + (3 \times 0) + (6 \times 0) \\ &= 2 \end{aligned} \begin{aligned} P(2) &= 2 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + 3 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + 6 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= (2 \times 0) + (3 \times 1) + (6 \times 0) \\ &= 3 \end{aligned} \begin{aligned} P(3) &= 2 \times \frac{(x - 2)(x - 3)}{(1 - 2)(1 - 3)} + 3 \times \frac{(x - 1)(x - 3)}{(2 - 1)(2 - 3)} + 6 \times \frac{(x - 1)(x - 2)}{(3 - 1)(3 - 2)} \\ &= (2 \times 0) + (3 \times 0) + (6 \times 1) \\ &= 6 \end{aligned}

## References

The following resources have been invaluable for me to learn the concepts discussed in this article.

You should definitely give them a read if you want to dive deeper into the topic.