3D Extrapolation in python (basically, scipy.griddata extended to extrapolate)

Question

I am using the griddata function in scipy to interpolate 3 and 4 dimensional data. It works like a champ, except that it returns a bunch of NaNs because some of the points I need are outside the range of the input data. Given that N-d data only works with the "linear" mode interpolation anyway, it should be a snap to have griddata do an extrapolation instead of just returning NaN. Has anyone done this or found a workaround? To clarify: I have unstructured data, so I can't use any of the functions that require a regular grid. Thanks! Alex

Would it help to fill the points outside the range with some constant value? In that case, you could just specify the fill_value — Dhara
Also, are you sure you want to extrapolate? sometimes, getting out NaNs and knowing you are going out of range is a much better choice. I have used Univariate splines from scipy, it silently extrapolates and the results can be quite "off" — Dhara
My situation is: I measure some values at a few points, and need to then calculate values at a bunch of other points through inter/extrapolation. So a constant value, or NaN really don't help. I know how temperamental the splines can be, so I was thinking linear would be a safe bet. I would like something that works on N-d data though. — Alex
This is actually a bit more nontrivial than it seems in these triangulation-based methods. For each point outside the tesselation, you need to choose which simplex you extrapolate from, and also have a fast algorithm to find this simplex. Also, if you want to avoid discontinuities in the extrapolation, some care is needed. After that is sorted out, however, the linear extrapolation is a simple matter. — pv.
I am just going to do a 2nd run as nearest and use that to fill in the gaps. thanks! — Alex

np8 np8 · Accepted Answer · 2020-12-04T21:31:14

One possibility to interpolate & extrapolate data with 3, 4 or actually any dimensions is with scipy.interpolate.Rbf

The get_data() function and plot_3d() function are attached to the end for convenience.

Example data

The example data looks like this (fourth dimension, w, is shown with a color). The data is irregularly spaced and not gridded.

x, y, z, w = get_data(N=200)
plot_3d(x, y, z, w)

Interpolation & extrapolation in 3d

First, let's setup new x and y coordinates. To make this more interesting, let's extrapolate to minus x and minus y direction. This forms the new x and y range of interest.

xs = np.linspace(-10, 20) # some extrapolation to negative numbers
ys = np.linspace(-10, 20) # some extrapolation to negative numbers
xnew, ynew = np.meshgrid(xs, ys)
xnew = xnew.flatten()
ynew = ynew.flatten()

Interpolation with scipy.interpolate.Rbf. Now,

from scipy.interpolate import Rbf
rbf3 = Rbf(x, y, z, function="multiquadric", smooth=5)
znew = rbf3(xnew, ynew)

plot_3d(xnew, ynew, znew)

There can be as many variables/dimensions as you want. The first arguments (x,y) are treated as coordinates for the nodes. The last argument before the function argument are the "values" which are to be interpolated (now: z).
The function argument can be used to control how to values are interpolated. This will affect the results so play with it with your data..
The smooth parameter can be used to smooth some noise from the data. If smooth is zero, the result is interpolating; it will go through all your data points. If it is positive value, the data is smoother. This will affect the results so play with it with your data..
Below are the results and the extrapolation is of course badly off. This is just to demonstrate that extrapolation is possible. You might want to fine tune function and smooth for desired results. Usually, data should not be extrapolated "much" (like in this example)

Adding fourth dimension

It is also possible to interpolate & extrapolate to the fourth dimension. Here is how:

rbf4 = Rbf(x, y, z, w, function="thin_plate", smooth=5)
wnew = rbf4(xnew, ynew, znew)
plot_3d(xnew, ynew, znew, wnew)

I created another Rbf instance for the fourth dimension and I used the znew calculated with the rbf3 (interpolation in 3d).
I changed the function to "thin_plate" since it seemed visually to perform better with this dataset.
Here is how the results looks like:

Appendix: `get_data` and `plot_3d`

For testing purposes:

import numpy as np

def get_data():
    np.random.seed(100)
    N = 200
    maxval = 20
    x = np.random.random(N) * maxval
    y = np.random.random(N) * maxval
    z = x ** 2 + np.sqrt(y) * y - y ** 3 + np.random.random(N) + 18 * y ** 2 * 2
    w = x ** 2 - np.log(y + (x * y) ** 2)
    return x, y, z, w


def plot_3d(x, y, z, w=None, show=True):
    import matplotlib.pyplot as plt
    from mpl_toolkits.mplot3d import axes3d

    fig = plt.figure(figsize=(10, 6))
    ax = axes3d.Axes3D(fig)
    ax.scatter3D(x, y, z, c=w if not w is None else "b")
    plt.show()