I am looking for an easy way (formula or algorithm) to create a polynomial function that gives the arbitrary preset output for the first values of x. For instance, the desired output can be $y = 1, 2, 3, 4, 5, 6, 100$ for $x = 1, 2, 3, 4, 5, 6, 7$. I can create one of the many functions by hand as:

$y = x + (x-1)(x-2)(x-3)(x-4)(x-5)(x-6)\cdot \frac{93}{720}.$

But for other desired outputs it is really painstaking. For example, I am working now with the desired output $y = 2, 10, 12, 16, 17, 18, 19, 200, 300$ when $x = 1, 2, 3, 4, 5, 6, 7, 8, 9$.


$x$ is always a natural positive number ($x > 0$). It doesn't matter what happens for other values of $x$ (greater than the ones given).
- I am looking for a easy way, as I don't know / I don't have / I have no experience with math related software.
- Also for a polynomial (or similar simple) function, as my mathematical background is quite limited, so I don't know any gamma, delta or other functions, series, etc.
- This quest is simply for recreation, to show friends how there are different solutions to the problems of the type: "Find the next term in this succession"


Polynomials in the standard form (as a sum of monomials) are actually the 'wrong' sort of 'encoding' to do interpolation. This post explains how to very efficiently interpolate any polynomial function, without ever computing the polynomial coefficients! In short, we can easily express any given polynomial sequence as a Newton series by computing successive forward differences of the given sequence, and not only can we use this to efficiently compute subsequent terms in the sequence, we can also quickly convert this representation back to standard polynomial form, in the same way the linked post demonstrates for the sum of the first $n$ perfect cubes.

For your example the iterated forward differences are easily computed by hand:

 1,  2,  3,  4,  5,  6,100,659,...
 1,  1,  1,  1,  1, 94,559,...
 0,  0,  0,  0, 93,465,...
 0,  0,  0, 93,372,...
 0,  0, 93,279,...
 0, 93,186,...
93, 93,...

Which also means that you can use this forward difference technique (which works for any polynomial sequence such as $1,4,9,16,25,\cdots$) to justify that you can extend any given sequence polynomially as a reasonable pattern.

That said, if too few terms of a sequence are given (which is often the case), it is actually easier to describe it as simply a list with no pattern, because the intended pattern would have a description that is longer than just listing the few terms given.

  • $\begingroup$ Thank you but I don't quite understand your answer, as I am not familiar with the forward differences technique, Newton series, nor how to get a polynomial from the finite differences. I've read the post you link, and Wikipedia's "Finite difference" article and the example they show, and have tried to use it in my two examples. $\endgroup$
    – user596191
    Sep 23 '18 at 15:45
  • $\begingroup$ cont. It seems that, as I have understood the procedure, it works for my first example (1,2,3,4,5,6, 100), but when I try it with my 2nd (2, 10, 12, 16, 17, 18, 19, 200, 300) where I get negative differences, and it doesn't work for me, surely because I'm not doing it right, as I don't fully understand the method. Any further help would be appreciated as it looks like a very interesting method. $\endgroup$
    – user596191
    Sep 23 '18 at 15:47
  • $\begingroup$ @Leonardo: As defined in the linked post, the forward difference operation on a sequence $f$ is the sequence of $(f(n+1)-f(n))$ over all indices $n$. Please study the example of cubes given there as well. If you still have a problem, show me your working for your second example. As long as you perform the forward differences until you reach a constant sequence, you can extend the original sequence just by extending that constant sequence and extending all the other sequences to obey the difference property, and the result will always be the least-degree polynomial that matches. $\endgroup$
    – user21820
    Sep 23 '18 at 15:52
  • 1
    $\begingroup$ @Blue: Yes that's right; Lagrange interpolation works for all finite sets of points that you want to interpolate with a polynomial. Newton series is just a 'magical' method that works for the special case of sequences. =) $\endgroup$
    – user21820
    Sep 23 '18 at 16:22
  • 1
    $\begingroup$ @user21820: I think I finally got it ! Thank you again ! $\endgroup$
    – user596191
    Sep 23 '18 at 18:07

As mentioned in a comment, the key is "Lagrange Interpolation". (See, for instance, Wikipedia or MathWorld.) The formal description can be a little opaque, but the idea is actually pretty straightforward: Build a polynomial out of parts, most of which go away with each input, then make the remaining part do what it needs to do.

For example, suppose you want the three (distinct) inputs $a$, $b$, $c$ to give the three outputs $p$, $q$, $r$. Lagrange claims that the polynomial you want has the form

$$h(x-b)(x-c) \;+\; j(x-c)(x-a) \;+\; k(x-a)(x-b) \tag{$\star$}$$

where $h$, $j$, $k$ are some as-yet-undetermined constants.

When $x=a$, two terms of $(\star)$ vanish, leaving $h(a-b)(a-c)$; likewise, $x=b$ yields $j(b-c)(b-a)$, while $x=c$ yields $k(c-a)(c-b)$. So, just choose $h$, $j$, $k$ to give the outputs you desire:

$$\begin{align} h(a-b)(a-c) = p \qquad\to\qquad h = \frac{p}{(a-b)(a-c)} \\[4pt] j(b-c)(b-a) = q \qquad\to\qquad \,j = \frac{q}{(b-c)(b-a)} \\[4pt] k(c-a)(c-b) = r \qquad\to\qquad k = \frac{r}{(c-a)(c-b)} \end{align}$$

Thus, we can re-write $(\star)$ informally as $\require{cancel}$ $$p\;\frac{\cancel{(x-a)}(x-b)(x-c)}{\cancel{(a-a)}(a-b)(a-c)} \;+\;q\;\frac{(x-a)\cancel{(x-b)}(x-c)}{(b-a)\cancel{(b-b)}(b-c)} \;+\;r\;\frac{(x-a)(x-b)\cancel{(x-c)}}{(c-a)(c-b)\cancel{(c-c)}} \tag{$\star\star$}$$

Clearly, this kind of thing can be extended to any number of inputs and outputs.

  • 1
    $\begingroup$ Sweet ! I was struggling trying to understand Lagrange Interpolation in Wikipedia, as suggested by @Gabriel Romon, and you've made it crystal clear. $\endgroup$
    – user596191
    Sep 23 '18 at 11:38
  • 2
    $\begingroup$ And a nice thing about this method is you perhaps don't need to find the actual coefficients of the final polynomial for your purposes. Write out the solution as a sum of $n$ fractions, maybe substituting in the $a$, $b$, $c$, ... variables, and say "we could multiply all that out and add up like terms, but it is clearly a polynomial and will give the specified sequence when $x=a$, $x=b$, etc.". $\endgroup$
    – aschepler
    Sep 24 '18 at 1:18

As Gabriel Romon commented: Lagrange interpolation does the trick.

The polynomials form a vector space (which works very similar to the 3D coordinate space) where the coefficients work as coordinates and the monomials including $x^0=1$ as orthogonal standard base. The Lagrange polynomials form another one. They are constructed from your given sampling points ($x$-values). The nice thing is that each Lagrange polynomial is zero in all but one sampling points. And in this sampling point it has the value $1$. So you can multiply each of the Lagrange polynomials with the function value at the sampling point where it has a nonzero value and sum over all these polynomials.

This looks as easy as $$ y = \sum_{i=1}^n y_i \cdot \ell_i(x). $$ As you can see there is no $x_i$ in this formula any more. This makes it possible to precompute the Lagrange polynomials $\ell_i$ and solve all your problems together in a secound step, as long as all functions are defined using the same sampling points.

  • 1
    $\begingroup$ That would save a lot of time! A good extension to @Blue answer! Great ! $\endgroup$
    – user596191
    Sep 23 '18 at 11:59

In a similar manner to user21820, you can take successive differences, so with your $y = 2, 10, 12, 16, 17, 18, 19, 200, 300$ example you would get something like this

2   10  12  16   17 18   19 200   300   
     8   2   4    1  1    1 181   100   
        -6   2   -3  0    0 180   -81   
             8   -5  3    0 180  -261   
                -13  8   -3 180  -441   
                    21  -11 183  -621   
                        -32 194  -804   
                            226  -998   

and if you assume that the final $-1224$ continues along the sequence you could extrapolate to get something like this:

2   10  12  16   17 18   19 200   300   -4030   -31346  -135672 -444797
     8   2   4    1  1    1 181   100   -4330   -27316  -104326 -309125
        -6   2   -3  0    0 180   -81   -4430   -22986   -77010 -204799
             8   -5  3    0 180  -261   -4349   -18556   -54024 -127789
                -13  8   -3 180  -441   -4088   -14207   -35468  -73765
                    21  -11 183  -621   -3647   -10119   -21261  -38297
                        -32 194  -804   -3026    -6472   -11142  -17036
                            226  -998   -2222    -3446    -4670   -5894
                                -1224   -1224    -1224    -1224   -1224

suggesting the next terms could be $-4030,-31346,-135672,-444797$ (that looks slightly implausible with the dramatic drop to negative values, and this uncontrolled behaviour is why polynomial extrapolation is generally seen as dangerous)

This also points towards finding the lowest degree polynomial fitting your data: there is an eighth degree polynomial where the coefficient of $x^8$ is $\frac{-1224}{8!} = \frac{-17}{560}$. So if you consider $y-\frac{-17}{560}x^8$ and take the differences again then the final row will be zero and the penultimate row $5734$ suggesting the coefficient of $x^7$ is $\frac{5734}{7!}=\frac{2867}{2520}$. Continue in this way with $y-\left(\frac{-17}{560}x^8 + \frac{2867}{2520}x^7\right)$ and so on, and you get the lowest degree polynomial fitting your data being

$y=-\frac{17 {{x}^{8}}}{560}+\frac{2867 {{x}^{7}}}{2520}-\frac{143 {{x}^{6}}}{8}+\frac{55169 {{x}^{5}}}{360}-\frac{187277 {{x}^{4}}}{240}+\frac{432067 {{x}^{3}}}{180}-\frac{362567 {{x}^{2}}}{84}+\frac{143421 x}{35}-1536$

  • 2
    $\begingroup$ Thanks for taking the time to actually compute that example for the asker. And your comment that polynomial extrapolation is dangerous is also spot on. Even interpolation is dangerous; roots can appear and disappear in the interpolated range with very slight perturbation of the points, and absolute value can be changed by an arbitrarily large amount if two points started sufficiently close. $\endgroup$
    – user21820
    Sep 23 '18 at 16:25
  • $\begingroup$ Thank you for your elaborated answer. I had gotten to a table like the first you have posted, and I can follow easily to the second and how to get the coefficient for x⁸. I am afraid, however, that I am lost trying to learn how to get the other coefficients $\endgroup$
    – user596191
    Sep 23 '18 at 17:01
  • 1
    $\begingroup$ @Leonardo Take $y-\frac{-17}{560}x^8$, i.e. $y+\frac{17}{560}x^8$ and then redo the difference table and see if you get the $5734$ I mentioned $\endgroup$
    – Henry
    Sep 23 '18 at 19:09

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy