Lagrange interpolation error formula

Step 1 (a clever auxiliary function). Fix a point $x$ that is not one of the nodes (if it is, the error is trivially $0$). Let $w(t) = \prod_{i=0}^{n} (t - x_i)$ and define the constant $c = \dfrac{f(x)-P(x)}{w(x)}$ (well-defined since $w(x)\ne0$). Define the auxiliary function $g(t) = f(t) - P(t) - c\, w(t)$.

Step 2 (count the roots of g). At every node $x_i$, both $f(x_i)-P(x_i)=0$ (by the interpolation property) and $w(x_i)=0$ (by definition of $w$), so $g(x_i)=0$ for all $n+1$ nodes. Also, by the choice of $c$, $g(x) = f(x)-P(x) - c\,w(x) = f(x)-P(x) - [f(x)-P(x)] = 0$. So $g$ has $n+2$ distinct roots: the $n+1$ nodes plus $x$ itself.

Step 3 (apply Rolle's theorem repeatedly). Between each pair of consecutive roots of $g$ (there are $n+1$ such gaps among $n+2$ roots), Rolle's theorem gives a point where $g'$ vanishes, so $g'$ has at least $n+1$ roots. Repeating this argument on $g', g'', \dots$ loses one root each time a derivative is taken, so after $n+1$ applications, $g^{(n+1)}$ has at least one root $\xi$ in the interval.

Step 4 (differentiate and solve for the error). Since $P$ has degree at most $n$, its $(n+1)$-th derivative is $0$; and $w$ is a monic polynomial of degree $n+1$, so $w^{(n+1)}(t) = (n+1)!$ for every $t$. Differentiating $g$ gives $g^{(n+1)}(t) = f^{(n+1)}(t) - 0 - c\,(n+1)!$, and setting $g^{(n+1)}(\xi)=0$ gives $c = \dfrac{f^{(n+1)}(\xi)}{(n+1)!}$. Recalling $c=\dfrac{f(x)-P(x)}{w(x)}$ and solving for the error gives exactly $f(x) - P(x) = \dfrac{f^{(n+1)}(\xi)}{(n+1)!} \prod_{i=0}^{n} (x - x_i)$.