Get polynomial coefficients from interpolation splines in R

1.3k Views Asked by At

I have a set of measured values that I'd liked to interpolate in R using cubic splines. Since these are just piecewise polynomials I'd subsequently like to integrate the interpolation function algebraically. Therefore I need the coefficients. Is there a way to obtain these?
Calling splines::interpSpline(foo, bar)$coef doesn't seem to return the actual polynomial coefficients.

1

There are 1 best solutions below

0
On BEST ANSWER

The output of splines::interpSpline(x,y)$coef gives polynomial coefficients of the part between x(i) and x(i+1) in terms of the powers of (x-x(i)), not of powers of x. This makes sense because the resulting coefficients are of reasonable size and are easier to interpret: for example, each constant term is just y(i), the quadratic coefficient gives the concavity at x(i), and so on.

For example, this output

> x <- c(1,3,6,9)
> y <- c(3,1,4,1)
> splines::interpSpline(x,y)$coef
     [,1]        [,2]       [,3]        [,4]
[1,]    3 -1.54054054  0.0000000  0.13513514
[2,]    1  0.08108108  0.8108108 -0.16816817
[3,]    4  0.40540541 -0.7027027  0.07807808
[4,]    1 -1.70270270  0.0000000  0.00000000

means that

  • on the interval [1,3] the polynomial is 3 - 1.54054054*(x-1) + 0.13513514*(x-1)^3
  • on the interval [3,6] the polynomial is 1 + 0.08108108*(x-3) + 0.8108108*(x-3)^2 - 0.16816817*(x-3)^3
  • on the interval [6,9] the polynomial is 4 + 0.40540541*(x-6) - 0.7027027*(x-6)^2 + 0.07807808*(x-6)^3

I don't see much use of the last line, which describes linear continuation of the spline beyond x=9, the right endpoint of the data.

To integrate these is no more difficult than to integrate powers of x, but of course one needs to choose the constants of integration if the goal is to get a continuous antiderivative. The choice of the form of polynomials makes dealing with constants of integration easier. Assuming we pick the antiderivative that has value 0 at the left endpoint, the rest is as follows:

  • on the interval [1,3] the antiderivative is 3*(x-1) - 1.54054054*(x-1)^2/2 + 0.13513514*(x-1)^4/4
  • on the interval [3,6] the antiderivative is C1 + 1*(x-3) + 0.08108108*(x-3)^2/2 + 0.8108108*(x-3)^3/3 - 0.16816817*(x-3)^4/4. Here C1 is the value of the previous antiderivative at x=3.
  • on the interval [6,9] the antiderivative is C2 + 4*(x-6) + 0.40540541*(x-6)^2/2 - 0.7027027*(x-6)^3/3 + 0.07807808*(x-6)^4/4. Here C2 is the value of the previous antiderivative at x=6.