Polynomial Interpolation with Matlab polyfit

Given $N$ data points, does polyfit of degree $N-1$ produces the unique interpolating polynomial?

For concreteness, here is a code example:

x=[1:10]
y= x.^3;
pp = polyfit(x,y,9)

where the following warning is issued:

distinct X values, reduce the degree of the polynomial, or try
centering and scaling as described in HELP POLYFIT.
> In polyfit (line 79)

pp is indeed $x^3$, as it should be, but why is this ill-conditioned? Interpolation is only ill-conditioned with respect to some bases (say, the monomials $1,x,x^2,\ldots$), but not with respect to others?

According with the documentation polyfit is a function that it can make an approximation in the sense of least square. For this reason you can choose the last input parameter, i.e. the degree. Oblivious if you put the degree equals to $N-1$ you are in the interpolation environment (sorry for the pedantic part :-) ).

In this case, as Federico's answer explain, the function build the Vandermonde matrix and solve the relative linear system (see here). This matrix is ill-conditioned and when matlab try to run

p = V\y

then here come the problem. Note that matlab try to use a $QR$ decomposition to solve the linear system and when the matrix $R$ is ill-conditioned with 1-norm matlab/octave, $||R||_1 ||R^{-1}||_1 > 10^{10}$ , there is the error message. I find this in my course notes in Italian, at the moment I can't get to you a better reference than this (see slide 10).

To avoid this issue you can use Barycentric Lagrange Interpolation, see this pdf by Berrut and Trefethen. The abstract is a good resume:

Barycentric interpolation is a variant of Lagrange polynomial interpolation that is fast and stable. It deserves to be known as the standard method of polynomial interpolation.

You can find the respective code, by Greg von Winckel, here.

• Thanks. Another option would be to use another basis, e.g., Legendre Jan 20 '18 at 16:30

In the end, this is a linear system. Its matrix is ill-conditioned: cond(vander(1:10)) returns 2.1063e+12. So if you change y by a certain amount, the change is amplified (in the worst case, using relative errors) by a factor 2.1063e+12.

The picture should be clear even without formulating statements such as "Interpolation is only ill-conditioned with respect to some bases" (which is not completely correct: if your nodes are the roots of unity, for instance, interpolation in the monomial basis is perfectly well conditioned).

In any case, here the interpolation basis is fixed: the function computes the coefficients of the interpolating polynomial with respect to the monomial basis, so that's what we have to return in the end.

• Thanks! Per your remark about the nodes, you're entirely correct, and I wasn't precise. However, if I want to compute the interpolating polynomial, is there a built-in option to do it in a well-posed basis? Jan 20 '18 at 15:26
• And per the original question, is it the case that polyfit gives the least squares polynomial, and in the special case of number of samples = degree of polynomials, it gives an interpolating polynomial? Jan 20 '18 at 15:27
• @AmirSagiv That's sort of a moot answer, but every interpolation problem is perfectly well conditioned in the Lagrange basis, because the solution is the input vector $y$ itself. Jan 20 '18 at 15:47
• @AmirSagiv Yes, in the case length(x)=degree+1, the least squares residual is 0 (assuming x has distinct entries). Jan 20 '18 at 15:48