Search code examples
pythonscipyinterpolationpythonxy

Get the formula of a interpolation function created by scipy


I have done some work in Python, but I'm new to scipy. I'm trying to use the methods from the interpolate library to come up with a function that will approximate a set of data.

I've looked up some examples to get started, and could get the sample code below working in Python(x,y):

import numpy as np
from scipy.interpolate import interp1d, Rbf
import pylab as P

# show the plot (empty for now)
P.clf()
P.show()

# generate random input data
original_data = np.linspace(0, 1, 10)

# random noise to be added to the data
noise = (np.random.random(10)*2 - 1) * 1e-1

# calculate f(x)=sin(2*PI*x)+noise
f_original_data = np.sin(2 * np.pi * original_data) + noise

# create interpolator
rbf_interp = Rbf(original_data, f_original_data, function='gaussian')

# Create new sample data (for input), calculate f(x) 
#using different interpolation methods
new_sample_data = np.linspace(0, 1, 50)
rbf_new_sample_data    = rbf_interp(new_sample_data)

# draw all results to compare
P.plot(original_data, f_original_data, 'o', ms=6, label='f_original_data')
P.plot(new_sample_data, rbf_new_sample_data, label='Rbf interp')
P.legend()

The plot is displayed as follows:

interpolation-plot

Now, is there any way to get a polynomial expression representing the interpolated function created by Rbf (i.e. the method created as rbf_interp)?

Or, if this is not possible with Rbf, any suggestions using a different interpolation method, another library, or even a different tool are also welcome.


Solution

  • The RBF uses whatever functions you ask, it is of course a global model, so yes there is a function result, but of course its true that you will probably not like it since it is a sum over many gaussians. You got:

     rbf.nodes   # the factors for each of the RBF (probably gaussians)
     rbf.xi      # the centers.
     rbf.epsilon # the width of the gaussian, but remember that the Norm plays a role too
    

    So with these things you can calculate the distances (with rbf.xi then pluggin the distances with the factors in rbf.nodes and rbf.epsilon into the gaussian (or whatever function you asked it to use). (You can check the python code of __call__ and _call_norm)

    So you get something like sum(rbf.nodes[i] * gaussian(rbf.epsilon, sqrt((rbf.xi - center)**2)) for i, center in enumerate(rbf.nodes)) to give some funny half code/formula, the RBFs function is written in the documentation, but you can also check the python code.